Multi-Zone Shell Model for Turbulent Wall Bounded Flows 
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We suggested a Multi-Zone Shell (MZS) model for wall-bounded flows accounting for the space 
inhomogeneity in a "piecewise approximation" , in which cross-section area of the flow, S, is subdi- 
vided into "J-zones". The area of the first zone, responsible for the core of the flow, Si — S/2, and 
areas of the next j-zones, Sj, decrease towards the wall like Sj oc 2~ J . In each j-zone the statistics 
of turbulence is assumed to be space homogeneous and is described by the set of "shell velocities" 
u n j(t) for turbulent fluctuations of the scale oc 2~ n . The MZS-model includes a new set of complex 
variables, Vj(t), j = 1, 2, . . . oo, describing the amplitudes of the near wall coherent structures of the 
scale Sj ~ 2~ J and responsible for the mean velocity profile. Suggested MZS-equations of motion for 
Unj(t) and Vj(t) preserve the actual conservations laws (energy, mechanical and angular momenta), 
respect the existing symmetries (including Galilean and scale invariance) and account for the type 
of the non-linearity in the Navier-Stokes equation, dimensional reasoning, etc. The MZS-model 
qualitatively describes important characteristics of the wall bounded turbulence, e.g., evolution of 
the mean velocity profile with increasing Reynolds number, He, from the laminar profile towards 
the universal logarithmic profile near the flat-plane boundary layer as Tie — > oo. 

PACS numbers: 05.45.-a, 47.27.-i, 47.27.Nz, 47.27.Eq, 47.60.+i 
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INTRODUCTION 
A. Background 

Three simple turbulent flows - in a channel, in a pipe 
and near a flat plane - play a prominent role in our 
understanding of spatially inhomogeneous wall bounded 
flows, similar to the celebrative role of the developed 
homogeneous turbulence in the understanding the uni- 
versal statistical behavior of fine-scale turbulence. On 
the long way of understanding homogeneous turbulence 
there appeared various phenomenological cascade models 
of turbulence (Richardson-Kolmogorov-41 concept of tur- 
bulence, Kolmogorov-62 log-normal and "multi-fractal" 
models of intermittency), many closure procedures (like 
Kraichnan Direct Interaction Approximation), various 
field theoretical approaches. Last but not least we have 
to mention so-called "shell models" of turbulence (like 
GOY shell model [1, 2] together with its "Sabra" im- 
provement [3] , and many others [4] ) . Separately we want 
to mention Zimin's shell model (V. Zimin, private com- 
munication, see also [5, 6]) which has been derived from 
the Navier-Stokes equation (NSE) using a vector wavelet 
decomposition [7] of the velocity field, and which involves 
no empirical or ad hoc parameters. 

The shell models are systems of ordinary differen- 
tial equations which mimic the statistically homogeneous 
isotropic turbulent velocity field in some interval of scales 
(say, within some "shell" in the Fourier space) by one or 
few "shell velocities" u n (t) [4], The shell models have 
the same (quadratic) type of nonlinearity, as the NSE, 
respect the conservation of energy (in unforced, inviscid 
limit) and have the build-in "locality" of interaction of 
neighboring scales, reflecting scale- by-scale energy trans- 
fer toward dissipative scales. Surprisingly, the shell mod- 
els allows to mimic almost everything we know (experi- 
mentally, theoretically or by direct numerical simulation) 



about highly nontrivial statistics of fine-scale turbulence. 
This include, for instance, the "intermittent" behavior of 
the velocity structure functions (which are simultaneous, 
two-point nth-order correlation function of velocity dif- 
ferences), the fusion rules (which govern asymptotic be- 
havior of many-point velocity correlation functions), and 
so on. A possible reason for such a lucky success is that 
the above mentioned (and some other) characteristics of 
the turbulent statistics are robust and depend only on 
some very general physical requirements, such as respect 
of the actual conservation laws , the scale invariance, the 
type of nonlinearity in the NSE, and the locality of inter- 
action. All these features are accounted for in the shell 
models. 

Unfortunately, the shell models in their traditional for- 
mulations describe only the space homogeneous turbu- 
lence, leaving aside the wall bounded turbulence, which 
play a much more important role in practical applica- 
tions. 

Turbulent flows at high Reynolds numbers, Tie 3> 1, 
contain such a wide range of excited lengths and time 
scales that the direct numerical simulations (DNS) of the 
NSE are impossible for the foreseeable future. Conse- 
quently, practical engineering calculations are based on 
some model simplifications of NSE, with Reynolds stress 
models being the most popular approach, see e.g. the 
book [8], review [9] and references therein. The idea of 
Osborne Reynolds (see, e.g. Ref. [10]) was to divide the 
velocity field into the mean flow part V(r) and the tur- 
bulent fluctuating part u(r,t) with zero mean and to 
approximate in some way the hierarchy of equations for 
various correlation functions (correlators). Equation for 
V(r) contains so-called Reynolds stress term, a second 
order correlator of u(r,t). The right hand side (RHS) 
of equation for the Reynolds stress contains five differ- 
ent terms: the rates of production, dissipation, turbulent 
transport and viscous diffusion, and the velocity pressure 
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gradient term. These are one-point, second and third or- 
der correlators of velocity and velocity gradients and the 
pressure-gradient correlator. The equation for only one 
such object, the dissipation rate, already contains eight 
correlators up to the 4th order in velocity, which are usu- 
ally modeled by various closure procedures in terms of 
the lower order objects. The simplest, old-fashion Mil- 
lionshchikov's closure (Monin and Yaglom, [10], p. 241) 
is often invoked. It uses a non-realistic assumption of 
Gaussian statistics of turbulence. To improve the sit- 
uation one can use a set of phenomenological constants, 
which can be found by comparison of the results of model 
calculations with the results of DNS or experiments in the 
benchmark flows. 

There where attempts to use more advanced field the- 
oretical approaches developed in the theory of homoge- 
neous turbulence, e.g. the Yakhot-Orszag version of the 
Renormalization Group (RNG) approach [11, 12]. In- 
stead of going into detailed criticism of this approach, 
already done in Ref. [13], we just make two general re- 
marks to the RNG approach, that are also relevant to the 
most attempts of a straightforward transfer of the field 
theoretical methods of fully developed, fine scale, homo- 
geneous turbulence to the case of wall bounded flows. 
First, in most cases, (RNG, diagrammatic perturbation 
approach, etc.) the turbulence is assumed to be excited 
by some artificial external force with Gaussian statistics. 
This is reasonable simplification of the real picture, if one 
deals with turbulent scales that are deep enough in the 
inertial interval. However, this is definitely not a real- 
istic assumption for large scales, that are important in 
the energy and mechanical momentum balance in wall 
bounded turbulence. Second, in most cases, the field- 
theoretical approaches to homogeneous turbulence are 
usually formulated in the fc-representation or explicitly 
assumed the space homogeneity. In this way one gets 
required closure relationships, say between the effective 
turbulent viscosity u T , the density of the kinetic energy 
£ and the rate of energy dissipation (used in the popu- 
lar K-e version of the Reynolds stress model). However, 
at least two of these objects (v T and £ ) are not locally 
defined, they are dominated by the largest eddies in the 
system, usually of scales close to the distance to the wall. 
Therefore one has to be extremely careful applying the 
resulting relations to the wall bounded flows in which the 
characteristic length of inhomogeneity is exactly the dis- 
tance to the wall. A price to pay for such a simplification 
is that the phenomenological constants may depend on 
the flow geometry or even on the position in the flow. 

Introducing enough adjustable parameters (sometimes 
geometry dependent), one can reach an engineering goal 
to model in computer some mean and turbulent char- 
acteristics of particular flows of practical importance. 
However, an important aspects of the basic physics of 
wall bounded turbulence remain unclear, being masked 
by numerous details, or even incorrectly reproduced. 

The main goal of this paper is to suggest a physically 
transparent and analytically analyzable model of wall 



bounded flows, oriented mainly on physical community. 
The model describes an interplay of two main physical 
phenomena in the wall flows: the energy cascade toward 
small scale (like in the developed homogeneous turbu- 
lence) and the cascade of the mechanical momentum to- 
ward the wall in the physical space. Our model is a gen- 
eralization of the shell model of homogeneous turbulence 
on the case of inhomogencous turbulence and accounts 
for a non-uniform profile of the mean velocity. Simpli- 
fying assumptions are made from the very beginning, at 
the level of basic, dynamical equations of motion. These 
equations involve only two parameters, responsible for 
the energy and mechanical momentum fluxes. These pa- 
rameters can be evaluated from the NSE, but currently 
they are chosen to reproduce the Von Karman constant 
k k and constant B in the universal log-profile of the mean 
velocity near flat plane. 

Our model is oriented toward the classical examples of 
the wall turbulence, like the channel and the pipe flows, 
the planar and the circular Couette flows. The physical 
description of the model and its equations of motion arc 
presented in the following subsection. 



B. Multi-zone Shell Model in brief 

For dealing with inhomogeneous wall turbulence we 
suggest in this paper a piecewise homogeneity approxi- 
mation in which the cross-section area of the flow, S, 
is subdivided into a set of "j -zones". The area of the 
first zone, responsible for the core of the flow, Si w S/2, 
and the areas of the next j-zones, Sj, decrease towards 
the wall like Sj oc 2~ J . In each j-zone, the statistics of 
turbulence is assumed to be space homogeneous and is 
described by its own shell model 

dUn ^ = - v n K 2 n u nj +J\f nj +A nj w n , (0.1) 

for the "shell velocities" u n j(t), which are responsible 
for turbulent fluctuations of the (dimensionless) scale 
s n ~ 2~™, referred for below as (nj)- eddies. Hereafter 
A n j is the Kronickcr symbol (= 1 for n = j and = oth- 
erwise). Equation (0.1) accounts for the viscous damping 
term with some effective viscosity v n ~ v , where vq is 
the kinematic viscosity of the fluid. The effective "shell" 
wave vector K n cx l/s n . The nonlinear term in Eq. (0.1), 
N n j, is given by Eq. (1.17) and describes the usual triad 
interaction of nearest shells, inside of a given j-zonc, 
Eq. (1.13), and some inter-zone interaction term of sim- 
ilar type. The production term w n is responsible for the 
energy flux from the mean flow to the turbulent subsys- 
tem and is given below by Eq. (2.39). 

Our goal is to describe the mean velocity profile 
(V(p,t)), in which p is two dimensional radius-vector 
in the cross-section of the flow, p_Lx, x is the stream- 
wise direction. To this end we introduce additional vari- 
ables Vj(t) with a prescribed space dependence &j(p), 
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uniquely determined by the flow geometry. The vari- 
ables Vj(t) can be understood as complex amplitudes 
of the near- wall coherent structures of the dimensionless 
scale sj, which is the same, as the scale of (nj)-eddies: 
Sj = h n for n = j. The functions &j(p) are chosen such 
that 3>j+i(p) ~ &j(2p) and form an orthonormal, (but 
incomplete) basis. We call it "PM-basis" , because (in 
spite of its incompleteness) it is chosen such as to rep- 
resent exactly the densities of the mechanical linear and 
angular momenta, V and M. , in terms of Vj (t) only: 

V = £ Sj Re[Vj(t)] , M = £ SjRj Im[V*(t)] ■ (°- 2 ) 



Here Rj is the characteristic distance of the j-zone from 
the centerline of the flow. The PM-basis also allows one 
to reconstruct the spatial dependence of the mean flow 
(and its time depending fluctuations) with a finite (but 
very good) accuracy 

V(p,t) « V PM (p,t) = Re[2^(*)*i(P)] -(0-3) 



For the PM -velocities Vj(t) wc suggest a simple mo- 
mentum equation 



dt 



Vp + w. 
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(0.4) 



which includes the viscous term, VjK^Vj, pressure gradi- 
ent, Vp > 0, and the Reynolds stress term, Wj, that ac- 
counts for the exchange of the mechanical momentum be- 
tween nearest zones in the flow and is given by Eq. (2.38). 

It is crucially important that suggested multi-zone 
shell model, Eqs. (0.1) and (0.4), preserves (in the un- 
forced, inviscid limit) all the relevant in the problem inte- 
grals of motion: energy, E, mechanical momenta, V and 
M. 



C. The plan of the paper 

Section I is devoted to the statistical description of the 
"turbulent part" of the Multi-Zone Shell model (MZS- 
model), Eq. (0.1). First, in Sec. IA we describe a way 
from the NSE to a standard shell model of homogeneous 
turbulence, that allows generalization for space inhomo- 
geneous turbulence. Next, in Sec. IB we formulate a 
piece-wise homogeneity approximation: the cross-section 
area of the flow, S, is subdivided into set of "j-zones" , in 
each of them the statistics of turbulence is assumed to be 
space homogeneous. This allows us to use the standard 
shell model for every j-zone and to describe the turbu- 
lence in the whole flow by the set of "shell velocities" 
u n j{t) for turbulent fluctuations of the scale oc 1/2™. 

In Sec. II we derive the dynamical equation of motion 
for the so-called PM- velocities, Vj, that allows one to 
reconstruct with good accuracy the mean velocity pro- 
file V(p) in the cross-section of the flow. In particular, 



in Sec. II B we introduce the PM-basis for a wide class 
of wall bounded flows that connects V(p) and Vj. In 
Sec. IIC we derive Eq. (0.4) for Vj and for all involved 
in Eq. (0.4) terms. 

In Section III we discuss the conservation laws 
(Sec. IIIB) and the symmetries (Sects. IIIC and HID) 
of the MZS model. For the convenience of the reader it 
begins with the overview of the resulting MZS equations, 
presenting them, in Sec. Ill A, in the dimensionless form 
(3.3), convenient for the further analysis. In particular, 
we derive in Sec. HID the version of the MZS model for 
the turbulent boundary layer near flat plane, Eqs. (3.19). 

Section IV presents a detailed analytical study of the 
MZS equations in a set of approximations, realistic at 
various values of Tie. In Sec. IV A we show that the MZS- 
model describes the stable laminar velocity profile for 
small Tie and its instability at some threshold value of 
Tie = Tie„. Next, in Sec. IV B we study the MZS-model 
in the approximation of the near- wall eddies, neglecting 
the turbulent energy cascade. This effect is accounted for 
in Sec. IV C in the approximation of turbulent viscosity. 
Section IV D is devoted to numerical analysis of the MZS- 
model. 

In Section V, we summarize our findings and suggest a 
possible generalization of the model for the description of 
turbulent flows laden with the long-chain polymeric addi- 
tives (in connection with the problem of drag reduction) 
or with heavy micro-particles, etc. 



1. Some notations and definitions 

• x & p - the streamwise direction & the two- 
dimensional radius-vector in the cross-section of a flow, 
p_Lx. In a channel p = (y, z), with y as the wall-normal 
and z as the spanwise directions. 

• S±, Pi & L - the cross-section area, the perimeter 
and the characteristic length of a cross-section: 



S± = J dp, L = S ± /P ± . 



(0.5) 



In a channel of width 2H, L = H; in a pipe of radius R, 
L = R/2. 

• Vp - the external pressure gradient, that is a positive 
constant 



dx 



(0.6) 



• t & U T - the characteristic time & the velocity in 
the flow: 



T=^LjV~p, U T =^LVp~. 



(0.7) 



The wall shear stress is U^. 

• vo & Tie - the kinematic viscosity & the friction 
Reynolds number 



Tie = U T L/v 



(Oi 
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• /' & /" - the real and imaginary parts of some 
complex object / (constant, variable, function, etc.): 



/' = Rc{/} , f" = Im{/} . 



(0.9) 



• n, j & p are dummy indices (natural numbers), 
reserved for the scale (or shell), zone & position indices; 
some objects can be related both to the shells and zones. 
They are used with both n or j indices. 

• Sj & <jj - the fraction of the cross-section area, 
occupied by the j-zonc, & the fraction of the cross-section 
area occupied by all zones towards the wall, starting from 
the zone j: 



Sj 

Si = 



S± 



i=3 



= 1 



• Lj = sj L - the width of the j-zone, L n - the char- 
acteristic scale in the n-shell 

• n n = 1/(2 L n ) - the wave number in the n-shell; Kj - 
the wave number of the j-zone 

• £, P, & M - the densities of the energy, of the 
mechanical linear & angular momenta (x projections) 

• , p ± , & - the total rates of pumping (with 
+ ) and dissipation (with ~) of the conserved quantities 
£,V,kM 

• ef , pf, & - the rates of pumping (with + ) 
and dissipation (with ~) in the j-zone (of £, P, & M) 

• Ej , pj, & nij - the fluxes of £, P, & M from the 
j- to the (j + l)-zone 

• The scalar product of complex vector functions, 
A(p) & B(p): 



(A,B) = / A 



'(P)-B(p)^- 



(0.11) 



• 0+ (p) & 0~ (p) the even and odd eigenfunctions of 
the two-dimensional Laplace operator in the cross-section 
of the flow with no-slip boundary conditions, (— p) = 
±0±(p) 

• *j(p) = *j(p) +i&j(p) ~ -PM-basis functions, 
Eq. (2.4) 

• ^pm(p) & Vj - PM- velocity in the coordinate and 
j-representations, related by Eq. (2.14) 

the velocity of turbulent fluctuation of the 



scale cx 2 ™ in the j-zonc [(nj)-eddies]; Uj = u 
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STATISTICAL MULTI-ZONE SHELL MODEL 
FOR TURBULENT FLUCTUATIONS 



From NSE to shell models of homogeneous 
turbulence 



In this Subsection we present a re-derivation of a stan- 
dard shell model of space homogeneous turbulence in a 
way that allows us to generalize it in Sec. I B for the case 
of space inhomogeneity. 



1. "Cell basis", wavelets and "(np)- eddies" 

Consider for simplicity an incompressible turbulent ve- 
locity u(r, t) in a periodic box of size L x L x L. Instead 
of the r- or k- representation we introduce here a 11 cell ba- 
sis" ^ np (r), which is quite similar to wavelet bases (for 
an easy-to-read, introductory text about the theory of 
wavelets see, e.g. [14]). Similar to the wavelet bases, the 
cell bases reflect both spatial scales of turbulent struc- 
tures (like in the fc-representation) and their position in 
the physical space (like r-representation) but accounts 
for the actual boundary conditions of the flow. The cell 
index, (n,p), consists of the scale index n = 1,2, ...oo 
and the position index p = (p x ,p y ,p z ): 

The scale index n — 1, 2, . . . oo defines the character- 
istic width L n (in all directions) of the function *„ p (r), 
which some "(np)-cell" of the size L n = L/2 n "occupies" 
in space. 

The position index p = (p x ,p y ,p z ) defines the position 
R np of the (np)-cell: 



pL n , P = (p x , Py, P Z ) , 



(1.1) 



p a = 1,2,... T. 



a = x,y,z . 



One can imagine that, for any given n, the set of (np)- 
cells with p £ C n fills the periodical L 3 -box. The turbu- 
lent velocity field u(r,t) is given by: 

oo 

M M) = EE Rc{[/„ p (t)*„ p (r)} , (1.2) 

n=l peC„ 

where the amplitude of the "cell expansion", U np (t) is 
the velocity difference across the separation L n in the 
(np)-cell. 

It is convenient to normalize the cell-basis as follows: 

LLL 

JJJ K P (r)-^n^{r) dX< ^ 3 dZ = 2v n A nn ,A pp ,, 



000 
LLL 



^np(r)-^n' P '(r) dX f 3 dZ = ■ (1.3a) 



000 



Here v n — 2~ 3n is the dimensionless part of the total 
volume per one mode in the nth cell. 

The equations (1.2) and (1.3) give the Parseval identity 
for the density of the turbulent energy in the form: 



LLL 



£ = 



\u(r,t)\ 2 dxdydz ^ 



000 



71=1 



£n(t) = !)n E 



(1.4b) 



pec„ 



This equation supports our interpretation of U np (t) as 
the velocity difference across the separation L n in the 
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(np)-cell. We will refer to these fluctuations as to a (np)- 
eddy". 

A particular choice of the cell functions ^ np (r) is not 
important for us here. Notice only, that for large n ba- 
sic cell functions become scale invariant and may be ob- 
tained by dilatations of one (or few) n-independent func- 
tion \I r oc (x). In this limit the cell basis becomes the 
wavelet one with * 00 (x) as so-called "7^- wavelet" [14]. 
An explicit example of a divergence-free 3 dimensional 
vector 1Z- wavelet function given in [5]. 

The cell functions ^ n p(r) form a complete orthonor- 
mal basis, therefore one can derive the exact equation of 
motion for U np (t) by the Galerkin projection of NSE: 



dt 



^n,pp'U np ,(t) 



(1.5) 



p' 



+ y y T°°,°,, pp , p ,,u°, p ,(t)u°, lp ,,(t), 

n'n" P ' P " a'(r" — zn 

where a, a' ,o" = ± are sign indices and we accept the 
convention U~ p = U np and = U* p . The explicit 
forms of the damping parameters T and amplitudes T 
depend on the basis, see, e.g. [5]. 



2. Basic assumptions for the standard shell models 

Consider shortly the physical simplifications that al- 
lows one to reduce the NSE to the shell model of homo- 
geneous turbulence. Unlike the similar discussion in [5], 
we made an accent on the assumption of space homogene- 
ity and the possibility to relax this assumption in order 
to generalize shell models for the space inhomogeneous 
case. 

The standard shell models of homogeneous turbulence 
follow from the exact Eq. (1.5) with the following simpli- 
fying assumption of a statistical nature: 



U np (t) 



0, IX 



np\ 



1 . 



(1.6a) 
(1.6b) 



H6r6 A n p are time independent, random amplitudes, and 
denotes averaging over yet unknown statistics of 
A npi which is generated by the NSE (1.5). The physical 
meaning of u n (t) is the mean square velocity fluctuations 
of all (np)-eddies : 



u np (tW) = (KW| 2 > 



(1.7) 



with the mean velocity equal to zero (U np {f)j = 0, even 

if the time averaged shell velocity (u n (t)) ^ 0. 

The physical arguments behind Eq. (1.6) may be based 
on the fact that in the homogeneous turbulence all ve- 
locities U np (t) with different p have the same statistics. 
Eq. (1.6) therefore neglects only the difference between 
the actual time realizations of np- velocities U np (t) of the 



same scale (the same scale index n), but occupying dif- 
ferent cells (different position index p) . The ensemble of 
the time realizations is replaced by the time-independent 
ensemble of A np . With the assumption (1.6), Eq. (1.5) 
yields: 



du n (t) 



dt 



In 



ao'o" 
^nn'n" 



(t)+N n , (1.8a) 

^ ^ n,pp' A np A np r , (1.8b) 
P > 

E S ™nC'«n'(*K"(*). (1-8C) 
n'n" .a' a" 

^ ] T^nn'n" ,pp'p" A np A n ' p ' A n " p " . (1.8d) 



The correlation functions of A np can be evaluated with 
some reasonable statistical assumptions (for instance, 
with some closure procedure, like Direct Interaction Ap- 
proximation, see e.g. [5]). 

Under assumption (1.6), Eq. (1.4b) gives a usual equa- 
tion for the density of the total energy of the nth 

scale: £„(i) = ||m„(£)| • Different shell models corre- 
spond to various further simplifications of the nonlinear 
term (1.8c). For example, in the Sabra shell model [3]: 



N„ 



= i[aK n+1 u n+1 u n+2 

-CK rl _lM„_ 2 W„-l] 

oc 2™, a + b + c = 



bn n u 



nO n _iU n +i 



. 



(1.9a) 
(1.9b) 



Notice, that the scale index n in the Eqs. (1.8) and (1.9) 
for shell models becomes the shell index. 



B. Piece-wise homogeneity approximation and 
multi-zone shell model for turbulent fluctuations 

The turbulent fluctuations in the wall bounded flows 
are not space homogeneous due to the spatial dependence 
of the mean velocity profile, which, in its turn, is also af- 
fected by the turbulent fluctuations. Due to the inhomo- 
geneity of turbulence, the shell model approach discussed 
in Sec. I A have to be revised, which is the subject of this 
subsection. For concreteness we discuss the planar geom- 
etry of the channel flow of width 2H in the cross-stream 
direction y, see Fig. 1. In further analysis we consider 
only the low half of the channel, < y < H, having 
in mind that the flow in the second half of the channel, 
H < y < 2 iJ, is statistically identical to the that in the 
first one. 

Clearly, the core of the flow (say, for H/2 < y < H) 
may be approximately viewed as homogeneous. Let us 
call this region the "1-zone" . The next region H /A < y < 
H/2, in which the mean velocity profile V(y) begin to 
decrease towards the wall, we call the "2-zone" . Notice, 
that the width of the 2-zone, H 2 = H/4 is a half of 
the 1-zone width, H\ = H/2. Therefore, approximately 
with the same accuracy we can consider the statistics of 
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FIG. 1: Geometry of the channel and the plane Couette flow 
between two parallel planes separated by 2H in the cross- 
stream direction y. In the simple channel flow the pressure 
gradient is applied in the "streamwise" direction x and the 
mean velocity V(y) is oriented also along x. In the plane 
Couette flow the lower wall (y — 0) is moving in z (span- 
wise) direction with some velocity Vo, while the upper wall 
(y = 2H) is moving in the opposite direction. In this case 
V(y) has 2-projection. For both flows (and their hybridiza- 
tion) the three-dimensional velocity fluctuations are space- 
homogeneous in the x — z plane. 

turbulence in the 2-zone as homogeneous, but different 
from that in the 1-zone. 

Similarly, one expects that the mean velocity difference 
across each next j-zone of the width Hj = HjV will be 
more or less the same. This is the motivation to define 
" j -zone " in the scale invariant manner, as 2~ J H < y < 
2~ii- 1 )H and to approximate the turbulence inside each 
such zone as homogeneous. 

The approximation of the piecewisc homogeneity al- 
lows one to use the shell-model reduction, Eq. (1.6), in- 
side of each j-zone, similarly to that in the whole space 
for homogeneous turbulence: 

U np (t) u nj (t)A npj , l<j<n, (1.10a) 
~~A — n Ta |2 — i 
u nj = u nn , j >n . (1.10b) 

Here Pj belong to the j-zone, in the sense that the {npj)- 
cells are inside the j-zone. We introduced in Eq. (1.10a) 
the velocity of (nj)-eddies, nth shell in the j-zone. Equa- 
tion (1.10b) reflects the fact that the near-wall (nj)- 
eddies with zone index j > n belong simultaneously to 
the j-zone with j = n, sec Fig. 2. Therefore in our model 
all u n j(t) with j > n are just u nn (t). 

With Eqs. (1.10) one gets from Eqs. (1.4b): 

= E§M*), M*) = ^f^, (i.ii) 

where £ n j is the energy density of the nth shell in the 
j-zone of the width Hj = H/2^ . 

Using Eqs. (1.10) and (1.5) one gets the equation of 
motion of the multi-zone shell model for turbulent fluc- 
tuations: 

= -7n«nj(*) + My + A nj Wj , (1.12) 



in which we have added by hand the production term uij , 
describing the energy pumping to the turbulent system. 
This term will be clarified in the following Section by 
Eq. (2.39). 

The nonlinear term M n j in Eq. (1.12) describes the 
total energy balance for the nth shell in the j-zone. There 
are two distinct geometries: j < n and j > n. In the 
latter case it is enough to describe the energy balance for 
j = n, since all u n j — u nn for j > n. 

All eddies with j < n are fully placed in the same zone 
(e.g., the eddies of shells n — 2, 3 and 4, that belong to 
the 1st zone, see Fig. 2) . Therefore in this case we can 
use for M n j if) a standard shell model expression for N ni 
in which u n j =>■ u n . In this paper we adopt the Sabra 
version of N n j term, generalizing Eq. (1.9): 

Mnjf) = N nj {t), for j<n, (1.13a) 
N nj = i [a K n+ iu n +1 jU n+ 2,j (1.13b) 
+b 

The energy balance of the near- wall eddies (n = j) is 
quite different. As one sees in Fig. 2 (on the example of 
the eddy in the first shell), the near- wall eddies partici- 
pate in the triad interactions of three types: 

- Triad 1 involves one near- wall eddy and two bulk ed- 
dies [(un — U21 — uzi) triad of the lst-zone in the above 
example]. For this interactions we will use Eq. (1.9) but 
with different parameters: a\, b\ and c\\ 

- Triad 2 involves two near-wall eddies and one bulk 
eddy [(ui2 — u 2 2 — "32 ) triad of the 2st-zone in Fig. 2]. 
Here we will use Eq. (1.9) with the parameters: a 2 , b 2 
and c 2 ; 

- Triad 3 involves three near-wall eddies [(^13—1*23—1*33) 
triad of the 3rd-zone in the above example] . Here we will 
use again Eq. (1.9) with the parameters: a^, 63 and C3. 

The relationships between four set of interaction pa- 
rameters [(a, b, c) and ia p ,b p ,c p ), p = 1,2,3] may be 
found from: i) the requirement of the conservation of 
energy and ii) the "correspondence principle" : 
For the space homogeneous case (u n j is independent of 
the zone index j) the multi-zone model must coincide 
with the usual shell model for homogeneous turbulence 
(in our case, with the Sabra model). 

The above two requirements give: 

01 = a/2, a 2 = a 3 = a/4 , (1-14) 
61 - 6, b 2 = b 3 = b/2, 

c = ci = c 2 = c 3 . 

These equations may be interpreted as follows: 

- in the Triads 1 only the half of the largest eddy be- 
longs to the same zone as two smaller ones. This gives 
ax — a/2. Two smaller eddies in the Triad 1 fully belong 
to their zone. This corresponds to 61 = b and ci = c; 

- in the Triads 2 only the quarter of the largest eddy 
belongs to the same zone as two smaller ones. Therefore 
a 2 = a/4. In this triad only the half of the middle eddy 
belongs to the same zone as the smallest one. This corre- 
sponds to b 2 = 6/2. The smallest eddy in the triad fully 
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FIG. 2: Zones, shells and triad interactions in the multi-zone shell model for the channel geometry. Regions of localization of 
(np)-eddies are shown schematically as ellipses with corresponding numbers inside. The near-wall eddies have n = j. They 
also occupy all j'-zones towards the wall, with j > n. 



belongs to their own zone, and we take c 2 = c; 
- in the Triads 3 only the quarter of the largest eddy, the 
half of the middle one and the full smallest eddy belong 
to the same set of zones (Zones 3,4,.. in our example). 
This corresponds to the relationships 03 = a/4, 63 = 6/2 
and c 3 = c. 

The total contribution of three types of the triad in- 
teractions to the nonlinearity of the near-wall eddies can 
be summarized as follows: 

Nnn = ^ [ 2N n,n + N n>n+ 1 + iV„, n+2 ] . (1.15) 

One can join two Eqs. (1.13) and (1.15) and write: 



n+l 



(1.16) 



+-(iV n>n+2 - N nn ) 

Obviously, in the homogeneous case the second term in 
the RHS of Eq. (1.16) vanishes and one recovers the usual 
shell model for homogeneous turbulence with N n j => N n . 

Notice that the explicit form of this equation reflects 
the fact that we have accounted only for the triad inter- 
actions involving the nearest shells (n — 1, n and n + l) 
and the particular form of the channel subdivision on the 



zones: Hj = 2~ J H. This subdivision of the cross-section 
area on zones is reasonable for the scale invariant case 
of turbulent boundary layer near flat plane. The physi- 
cally motivated subdivision in particular flow geometries 
will be discussed in the following section. In general case 
Eq. (1.15) is changed as follows: 



i— n+l 



Mnj = N nj + A nj J2 ^(N m -N nn ), (1.17) 

oo 



with an arbitrary dependence of the (dimensionless) j- 
zone areas, Sj, on j and an arbitrary form the nonlinear 
term N n j. The fact that under the sum one has the dif- 
ference (N n j — N nn ) guarantees the correspondence prin- 
ciple, while the weights / <7j follow from the requirement 
of the conservation of the total energy, which in the case 
of arbitrary zone division is given by a natural general- 
ization of (1.11) with Hj/H replaced by sj; 



^n(^) — ^ SjSfij 



'■nj 



(1.18) 
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II. DYNAMICAL MULTI-SCALE MODEL FOR 
THE "PM- VELOCITY" 

A. Mechanical momenta, P and M, 
and "statistics vs. dynamics" dilemma 



makes it hardly possible to remain on the level of statis- 
tical description and one has to deal with the detailed, 
dynamical description of the PM -velocity V PM (p, t) in 
terms of the NSE. This is the subject of the following 
Subsections. 



The unbounded turbulence may be described in the 
reference system with a zero mean velocity. In that sys- 
tem the total linear momentum vanishes: P = 0. Due to 
the Galilean invariance the space homogeneous velocity 
does not interact with the turbulent fluctuations. There- 
fore one can consider turbulence in any other reference 
system with P^O with the same result for the statistic 
of turbulence. Hence, the mechanical momentum is not a 
relevant integral of motion for the unbounded turbulence. 
This is not the case for the wall bounded turbulence, in 
which the Galilean invariance is broken by the presence of 
the walls. The conservation laws for the total linear mo- 
mentum, P, (as well as for the total angular momentum, 
M) give an important constraint on the expected behav- 
ior of the system. For example, in the channel and pipe 
flows, the total input of the linear momentum due to the 
pressure gradient, acting in the whole cross-section area 
of the flow, must be equal to the dissipation of the mo- 
mentum on the walls due to the viscous friction. Clearly 
an analytical description of the wall bounded turbulence 
must respect the conservation of P and M (in the un- 
forced, inviscid limit). 

Consider the mean velocity field V(p) in the simple 
flows with translation symmetry in the streamwise direc- 
tion x||P, in which also M || x. Here p is the radius- 
vector in the cross-section of the flow, p_Lx. The exam- 
ples of such flows are the channel and the pipe flows, the 
planar Couctte flow, Fig. 1, the circular Couette flow, 
etc. Instead of P and M, it is more convenient to deal 
with their volume densities V and M, defined as follows: 



V = V x = (x, V) , 
M = M x = (K, V) , 



71 = x x p 



(2.1a) 
(2.1b) 



where (A,B) is the scalar product (0.11). 

Our idea is to divide the full velocity field V(p) into 
two parts, denoted as V PM {p) and V T (p), such that the 
"turbulent part" V T (p) does not contribute to V and M: 



(x,V T ) = 0, (K,V T ) = 0. 



(2.2) 



and to take a special care only on V PM {p) part, con- 
tributing to V and M.. This division may be done in 
many ways, our particular choice will be clarified below 
by Eq. (2.14). 

The description of V PM (p, t) may be statistical or dy- 
namical. The statistical description is the straightfor- 
ward "cell-expansion approach to shell models" , based on 
the NSE for V PM (p, t) in the cell representation, similar 
to Eq. (1.5). To continue, one has to find some reasonable 
statistical simplifications, similar to Eq. (1.10). However, 
since of V ^ and/or M. ^ 0, now also {V PM (p, t)) ^ 
and therefore A np cannot be approximated as zero. This 



B. PM-basis, properties and interpretation 

1. Construction of the PM-basis 

To construct a "physically motivated" basis for de- 
scription of the PM-velocities, consider two subsets 
of cigenfunctions of the Laplace operator, <f>m(p) an d 
</>~(p), satisfying the incomprcssibility condition. The 
no-slip boundary conditions are assumed in the cross- 
section of the flow, with the constraints 



Pm = (X,0+) ^ 0, 

TZ m = Ofc,0-)^O, 



(2.3a) 
(2.3b) 



dictated by Eq. (2.2). Notice that both sets 0+ and 0~ 
are chosen here as real and orthonormal. 
Introduce the complex PM-basis as follows: 



E 



(2.4a) 



m e Sj O 2 j ~ 1 < m < 2 3 ' - 1 . 



(2.4b) 



The normalization "radius" i? 7 -, is chosen such that 



(2.5) 



where 



*;.(p)^Rc[*,(p)], ^(pjElm^-lp)] 
Denote 

We show below, that 

oo 

E S J = 1 ' 



(2.6) 



(2.7) 



that allows us to understand Sj as a portion of the cross- 
section area, occupied by the j-zone. Then, the first of 
Eqs. (0.10) gives the definition of Sj, the area of j-zone, 
that is consistent with the requirement X^=i &j = ^-L- 
To find Rj from Eq. (2.5) compute also 



(&! = — V TZ 2 

3 m£Sj 



Then 



7 1 E K 



(2.8) 



(2.9) 
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2. Orthogonality and normalization conditions 

By construction (2.4-2.5) the PM-basis is orthogonal 
in the sense: 



(*J>*j') = 2s 3 A 33' > 



(2.10a) 
(2.10b) 



The idea behind the choice (2.4) is that the functions <I>^ 
and form the exact expansions of the uniform profile 
of unit height (<&'j) and of the linear profile (3>") : 



E*» = x < 



(2.11a) 
(2.11b) 



These relations can be easily proven by projection of both 
sides of Eq. (2.11) onto the Laplace basis 0^, using (2.3) 
and the definition of the PM-functions *j(p), Eq. (2.4). 

Let us show that Eq. (2.7) follows from the obvi- 
ous constrain (x, x) = 1 , after substitution x from 
Eq. (2.11a). Using also the orthogonality condi- 
tions (2.10) one gets 

oo 

i = (x )X )= j2 (*;•,*;•') 

^ oo oo 

= 7 E(*J+*J'*i'+*i') = E-J" 



We proved that all Sj add up to unity according to 
Eq. (2.7). 



3. PM -expansions 
Introduce the PM-projector 

P PM {f(p)} = 2ReJ2j§j^-f 3 iP) (2-12) 

= Re^T 1 ^,/)*^) , 

i 

and define the PM-velocity V PM (p,t) as the projection 
of the full field V(p,t) on the PM-basis: 

V PM (p)=P PM {V(p)} . (2.13) 
Using the normalization (2.10) one gets 

oo 

V PM (p,t) = Re{]TvS(t)$,(p)}, (2.14a) 

^•(t) = = (*,-, Vp M ) . (2.14b) 

The expansion coefficients Vj(t) can be understood as 
PMvelocities in the zone representation. They play the 



major role in our theory and will be called hereafter "the 
PM- velocities" or "zone velocities" . 

Substituting x and TZ from Eqs. (2.11) into the def- 
initions (2.1) of V and A4, one gets with the help of 
Eqs. (2.14) the PM-expansion of the linear in V inte- 
grals of motion: 



v = E(*j> v ) = E^'(*) : 



(2.15a) 



M = '£R j (*'>,V) = '£s j (R j Vj'(t)) .(2.15b) 
i=i j 

This is a proof of the statement that only PM part of 
the full velocity V(p,t), Eq. (2.14), contribute to the 
mechanical momenta. The turbulent part V T (p) does 
not contribute to the linear integrals of motion and will 
be considered as a part of the "turbulent ensemble" , de- 
scribed by the shell variables u n j(t). 

The energy density associated with the PM-velocity is 
as follows: 



£™(t) 



s± 



\V PM (P,t)\ 2 dp 
2 S± 



3 = 1 



(2.16) 



4- PM-basis for the channel and pipe flows 

The particular form of the PM-functions is geometry 
dependent. To give a feeling how these functions may 
look like, we discuss here two important cases, the chan- 
nel and the pipe flows. These two examples will serve us 
in the rest of the paper. Some properties of these specific 
bases are however more general and will be used in the 
derivation and analysis of the momentum equation for 
Vj(t). 

In the planar geometry, Fig. 1, p =>■ y and the functions 
4>m in Eq. (2.4) are given by 



<t>Uv) 
<t>™(y) 
4> m {y) 



Z02m(j/) , 

y/2 sm(k m y) , k 

2V2 



(2.17a) 



7r(2m - 1) 



nm/{2H), 
0.9/(2m- 1) . (2.17b) 



For the pipe geometry, 



0m (P) = x<Pm{p) , <t> m {p) = e<f> m (p) , (2.18a) 

<Pm{p) = MkrnP) / Jl(k m R) , (2.18b) 

J (k m R) = 0, Pm = 2/(k m R), (2.18c) 

where e is the polar- angle unit vector (e _L x, e _L p), 
Jo(0 and Ji(£) are the Bessel functions of the 0-th and 
1-st order, and R is the radius of the pipe. First of 
Eq. (2.18c) defines k m = £ m /P via zeros of the Bessel 
functions £ m : Jo (Cm) = 0. 
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FIG. 3: First functions <frj(£H) describing the mean flow with non-zero mechanical momenta in the channel of the width 2H. 
In the whole channel < £ < 2. Panels a and b: the PM-functions &j X (£H) (Panel a) and <&'j z (£H) (Panel 6) with j = 1, 2, 3, 4 
for the whole channel. Panels c and d: the PM-functions with j = 5, 6, 7, 8 in the near-wall region. Panel c: &j X (^H). Panel d: 



In Fig. 3 we plot &j^ x (y) and &^ z {y) for the channel as 
it follows from Eq. (2A) with Eq.'(2.17). The functions 

x (y) are symmetric with respect to the centerline of 
the channel y — H, whereas &j Z {y) are antisymmetric. 
For large j, &j x {y) and $" jZ (y) coincide in the lower 
half of the channel y < H and have an opposite sign for 
y > H . As it is clear from Panels c and d, in the near- 
wall (y <C H) region x (y) and z (y) almost coincide 
already for j = 5. 

For the channel flow one can find an explicit asymp- 



totic expression for <6j x (y) and 



oo (actually, j > 4). For y < H: 



in the limit 



~ $^(y)~$ un (2%V2ff), (2.19) 



*un(£) 



Si 



24 



Si 



(2.20) 



where Si(x) is the sine integral function. 



It is expected, that the asymptotic form of basic func- 
tions, $ M , is the same for any flow geometry, if one ex- 
press them as functions of the distance from the wall. 
For example, in the pipe one obtains (2.19), where y is 
the distance from the wall, y = R — p, and H is the ra- 
dius of the pipe, H = R. In Fig. 4 we show the collapse 
of the re-scaled functions <§>j(2H£/n 2 2 : >) for j = 5, . . . 8 
with the universal function $ U n(£) for the channel (Panel 
a) and pipe (Panel b, H = R) flows. 



5. Geometry of the j-zones in the channel and pipe flows 

a. Characteristic length of a flow. The cross-section 
of the flow can be characterized by two global parame- 
ters, the cross-section area, S±, and the length of its 
perimeter P±. One can organize from these two objects 
many combinations with the dimensions of the length, 
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FIG. 4: Asymptotic universality of the PM-basis functions. Panel a: Collapse of the re-scaled real part of the basis functions 
$ J (2_ff^/7r 2 2-') for n = 5.. .8 with the universal asymptotic function $ un (£) for the channel flow. The lines (from right to left) 
correspond to those in Fig. 3c. The left-most solid line denote <& un (£). Panel b: the same for the pipe flow, H — R. 



L((3) = P±(S±/Pj_)P , with arbitrary (3. For discussion 
of a turbulent flows under the external pressure gradient, 
the particular choice, (3 = 1, 



L = L(l) = S ± /P ± 



(2.21) 



is physically important. The reason is that the total ex- 
ternal accelerating force, applied on the unite length of 
the fluid (in the streamwise direction), is proportional to 
S±: 

F ac = VpS ± , 

while the total friction force is proportional to P± (under 
simplifying assumption of homogeneity along the perime- 
ter): 



F iri 



vqP. 



dV x (y) 



dy 



v=o 



The stationarity condition, F^ c — Ff r i c , allows one to re- 
late [dV x (y) /dy] y= o to the characteristic length, L defined 
by Eq. (2.21): 



dV x (y) 



dy 



y=o 



Vp 
i/ L 



In terms of the friction velocity, U T , Eq. (0.7), and the 
friction Reynolds number, Tie, Eq. (0.8), this gives the 
famous constraint for the wall bounded flows: 



dV x (y) 



dy 



= Ve T 



(2.22) 



In all mentioned global relationships the characteristic 
length L plays an important role. 

Notice that for the channel of the width 2H the length 
L is the distance from the centerline to a wall, L = H. 
For the pipe of the radius R the length L = R/2, which 
is twice smaller than the distance from the center to the 
wall. 



TABLE I: Parameters of the PM-basis for the planar and pipe 
geometries. 



3 


Planar 


Pipe 




Sj v 


Rj/H 


Vj/fo 


s j 


Sj v 


Rj/R 




1 


0.811 


1.62 


0.50 


6.48 


0.692 


1.38 


0.51 


2.77 


2 


0.122 


0.49 


0.77 


1.96 


0.185 


0.74 


1.03 


1.48 


3 


0.038 


0.30 


0.90 


1.22 


0.068 


0.54 


1.00 


1.09 


4 


0.015 


0.25 


0.95 


0.99 


0.029 


0.47 


1.00 


0.94 


5 


0.007 


0.22 


0.98 


0.89 


0.014 


0.44 


1.00 


0.87 


6 


0.003 


0.21 


0.99 


0.85 


0.007 


0.42 


1.00 


0.84 


7 


0.002 


0.21 


0.99 


0.83 


0.003 


0.41 


1.00 


0.82 


8 


0.001 


0.21 


1.00 


0.82 


0.002 


0.41 


1.00 


0.82 


oo 




2/tt 2 


1 


8/tt 2 




4/tt 2 


1 


8/tt 2 



b. Zones in the channel and pipe. In Tab. I we 
present parameters Sj for the channel and pipe geome- 
tries, that are given by Eq. (2.6) viap m . The parameters 
p m are taken from Eq. (2.17b) for the channel, and from 
Eq. (2.18c) for pipe geometry. 

Using Sj and L we introduce the characteristic width 
of the j -zone: 



(2.23) 



The idea behind this definition is that for the narrow 
near- wall zones Lj is exactly their geometrical width, Aj , 
i.e. the distance between the j-zone boundaries. 

Assuming that Aj is much smaller than the local curva- 
ture of the boundary, we evaluate the area of very narrow 
j-zone as Sj = Aj P± . On the other hand, this area is 
SjS±. Therefore, 



sjSJP ± 



for j ^> 1 
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TABLE II: Geometrical and characteristic widths Aj and Lj 



of jf-zone in the pipe geometry 



3 


1 


2 


3 


4 


5 


6 


Aj/R 


0.83 


0.105 


0.0356 


0.0148 


0.0071 


0.0035 




0.35 


0.093 


0.0340 


0.0145 


0.0070 


0.0035 



by projecting NSE (2.26) on the PM-basis, using 
Eqs. (2.10) and (2.14b) in the form: 



v j {t) = 2(* j ,v„ 



(2.28) 



1. The damping term 



Obviously, for the planar geometry, Aj = Lj for any j. 
To find the value of j, at which the zones become flat, 



consider the pipe geometry with Aj 



-i, where r,- 



is the radius of the circle, occupied by the first j-zoncs: 



1/2 



Pipe 



(2.24) 



Table II compares Aj/R, given by this equation with 
Lj/R = Sj/2. Clearly, for j > 3 these two parameters 
practically coincide. 



6. Interpretation of the PM -velocities 

The PM-expansions of the linear profiles (2.11) gives 
a simple interpretation of the velocities Vj (t) and Vj' (t) 
as the x- and e-projections of V PM (p, t) at some position 
within the j-zone, pf 

Vj(t) ^ V x ( Pj ,t) , V x (p, t) = x • V PM (p, t) , (2.25a) 
Vj'(t) V e ( Pj ,t) , V e (p, i) = e ■ y pM (p, t) , (2.25b) 

e = x x p/p . (2.25c) 

This point will be illustrated below in Fig. 14 by direct 
comparison of Vj and V(p) in appropriate coordinates. 

C. Simple momentum equation for Vj 



Consider the NSE for V PM (p, t) in the form: 
dV PM (p,t) 



dt 



-v AV PM +V P Z 



(2.26) 



-(V PM -V)V PM -(u-V)u. 



Here v§ is the kinematic viscosity, Vp = —dp/dx > 
is the pressure gradient in the streamwise direction x. 
The nonlinear term (V PM ■ V)Vp M describes the self- 
nonlincarity in the PM-subsystem and the last term, 
(it • V)w, is responsible for the effect of the Reynolds 
stress on the PM-velocity. Here 7TT is the ensemble av- 
eraging in the sense of Eq. (1.10). The nonlinear cross 
term, oc Vu (which has zero mean) is neglected. 

The goal is to get the equation of motion for the com- 
plex PM- velocities Vj(t): 



Projecting the NSE (2.26) on the real and imaginary 
parts of the PM-basis Eq. (2.4), and accounting only for 
the viscous term, oc i/q, one gets 



dV> 



-T'V' 



dV» 



-T"V! 

3 3 



dt 3 3 ' dt 

with, generally speaking, different 1^ and T": 



r' 

3 



-Vq- 



(2.29) 



However, in the pipe flow they are equal, 1^ = T" , since 
the functions and </>~ are the same, see Eq. (2.18). In 
the channel flow these functions, Eq. (2.17), are different 
and hence 1^ ^ T". Nevertheless, for large j, the basis 
functions ^ — > $" and Tj — » Y". Therefore it would 
be a reasonable simplification to neglect the (possible) 
difference between 1^ and T" and to write in Eq. (2.27): 



(damping) ■ = -TjVj 



(2.30) 



with the same damping term Yj, which is in between T^- 
and Y". For example, in the channel flow, Vj is respon- 
sible for the mean velocity profile V x (y) and thus "more 
important" than Vj' . In that case we will take Yj = Yj. 
In the planar Couette flow, the mean velocity profile is 
given by V z (y), which is connected to Vj', and thus one 
better take Y 3 = Yj. 

It is customary to represent the damping term via an 
effective viscosity, Vj, and an effective wave- vector, Kj, 
defined via the characteristic width of the j-zone Lj : 



2 Li 



Lj = Sj L . 



(2.31) 
(2.32) 



Parameters L and Lj were defined and discussed in 
Sec. IIB5b. The values of Sj and Vj for the channel 
and pipe flows are presented in Tab. I. 

Notice, that the damping term in Eq. (2.30) is diago- 
nal in j. This is a consequence of our definition of the 
PM-functions 4>j, Eq. (2.4): the functions with different 
j originate from different, and thus orthogonal eigenfunc- 
tions of the Laplace operator. 



2. The pressure term 



dV, 
dt 



= (damping)^ + (pressure) 



(2.27) 



+ (self- interaction) . 



( Reynolds-stress) ■ 



The equations (2.26) and (2.28) dictate: 

(pressure) . = 2Vf> (*,, x)/(*j, *j) . 



13 



Substituting x from Eq. (2.11a) into above expression, 
one gets for any flow geometry the same simple expres- 
sion: 



(pressure) • = Vp . 



(2.33) 



This means, that the uniform in p pressure gradient 
equally acts on all components V-{t). This remark- 
able result could be also obtained from the interpreta- 
tion of Vj, Sec. II B 6. Indeed, if we leave in the RHS 
of Eq. (2.26) only the pressure term, we get the p- 
independent V x (p,t) = Vpt, i.e. a homogeneous velocity 
profile. This means that Vj(t) is independent of j. 

3. The self-interaction term 



It is easily justified by noting that for j > 3 the j-zone 
contributions to M and M arc almost the same, see Ta- 
ble I for Rj. In the following we will drop the tilde in 
M., and write quasi- angular momentum simply as M.. 
Having in mind the conservation of the "complex mo- 
mentum" : 

C = V + iM = ^2sjVj, 



we suggest the following form of the Reynolds-stress 
term, Wj,: 

(Reynolds-stress) ■ = Wj = — (pj-i - pj) , (2.36) 



In various simple flows the self-interaction term in 
Eq. (2.26) is identically equal to zero due to geometri- 
cal constraint. It is so for the channel and the planar 
Couette flows, where 

(V ■ V)V => (V x A + V z -p\ V(y) = . 
\ ox oz J 

Here we skipped for the brevity the subscript " PM " . The 
same is true for the pipe geometry, where 

In the present paper we consider only flows with zero 
self-interaction: 



(self-ineraction)^ = . 



4- Effect of the Reynolds stress 



(2.34) 



In the PM-representation the Reynolds-stress term in 
Eq. (2.27) should be quadratic function of the turbulent 
velocities u n j and should scale as Kj = l/(2Lj) cx 2 3 . 
This term describes two physical processes, namely, (i) 
the exchange of the linear and angular momenta between 
different zones and, (ii) together with its counterpart Wj 
in the equation for u n j, the energy exchange between Vj- 
and u„j-subsystems. 

In the spirit of the shell models of turbulence, we con- 
sider here the form the Reynolds-stress term that ac- 
counts for the momenta exchange only between nearest 
j-zones ( j and j± 1) and preserves the relevant integrals 
of motions, the energy and the momenta. Notice, that 
in the turbulent channel and pipe flows with (M.) = 0, 
the conservation of V is much more important that the 
conservation of M. Therefore we can simplify the pos- 
sible structure of the inter-system interaction terms, Wj 
and Wj by using instead of Eq. (2.15b) for M its simpli- 
fied version, in which all Rj are assumed to be the same, 
Rj — » -Roo = const, and thus can be omitted: 



Quasi-angular momentum. (2.35) 



where pj is the momentum flux from the j- to the (j + 1)- 
zone. Indeed, this form provides the conservation of C in 
the inviscid, unforced limit: 

^ = s i W i = X^Pj" 1 - Pi) = • 

3 3 

For the momentum flux pj we also select the simplest 
form, assuming that only the near-wall turbulent eddies 
of the same scale ua = Uj (i.e. u n j with n = j) give a 



l 33 

nonzero momentum transfer: 

d 9 
^ - 2L U > 



(2.37) 



Here the dimensionlcss parameter d characterize the 
strength of the interactions. 

The equations (2.36) and (2.37) give us the following 
simple expression for Wj : 



Wj = -u?) • 



(2.38) 



The production term Wj , describing the energy pump- 
ing into j-zone of the turbulent subsystem, Eq. (1.12), 
is then uniquely determined by the requirement of the 
conservation of the total energy of the system: 



W j 



2a jL 



(Vj-V j+1 )u* . 



(2.39) 



5. Resulting momentum equation 



Collecting Eqs. (2.27), (2.30) 
finally gets: 

dVj(t)_ , 



(2.34) and (2.36) one 



dt 



= -UjK Z jV n + Vp + Wj 



(2.40) 



where Vj is the effective viscosity, given for the channel 
and the pipe flows in Table I, Kj is given by Eq. (2.32) 
and Wj - by Eq. (2.38). This equation together with 
Eq. (1.12) will be analytically and numerically analyzed 
in the following Sees. Ill and IV. 
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III. CONSERVATION LAWS AND 
SYMMETRIES IN THE MZS-MODEL 

A. Resulting dimensionless MZS-equations 

For the convenience of the reader we present here the 
full set of the MZS-equations which will be the subject of 
analytical and numerically studies. As the first step, we 
non-dimensionalize the MZS-equations, expressing time 
t and velocities Vj, u n j in the units t and U T : 



Here 





t 

1 


(3.1a) 




T 




Vj(rt) 


(3.1b) 


u nj (t) = 


u nj (rt) 

U T ' 


(3.1c) 


Vp = 


1, 


(3.1d) 



where time, r, and velocity, U T , are constructed from the 
"outer" characteristics of the flow, Vp and L: 

T = y/T/Wp, U T = ^/Vp~L. (3.2) 

In the rest of the paper we omit the tildes in the di- 
mensionless variables. The dimensionless MZS-equations 
take on the form: 

^ = -1^ + Vp 3 + Wj , V Pj = Vp = l(3.3a) 



du. 



dt 



(3.3b) 



j = 1,2,... oo , n = j,j + l,...oo. 



Here the shell variables the velocities of statisti- 

cally identical (nj)-eddies of characteristic scale s n that 
belong to j-zone of width Sj. Clearly, in Eqs. (3.3b) 
n > j. The PM- variables Vj(t) describe the velocities of 
coherent near-wall structures in the j-zone. In our no- 
tations the near wall turbulent (rtj)-eddies have n = j. 
These eddies occupy also all j-zones with j > n. There- 
fore in our approach 

u„ 3 = u nn = u n , for j >n . 

The terms in the RHS of the dimensionless MZS- 
equations (3.3) are as follows: 





Gj 2 9n 2 


1 

= 2V 


(3.4a) 




= dKj (uj_ 1 — Uj) , 




(3.4b) 


Wj 


= 2^(^' - V 3 + ^) U *3 ' °3 E 


= E s i' > 

?>3 


(3.4c) 


M n j 


= N n j + A nj ^ (N nj > - 
j'>3 j 




(3.4d) 














(3.4e) 




-\-b K n ll n _i jlln+ij — CK n — 


\U n -2,j u n- 





a + b + c = , ^ Sj = 1 



The equations (3.4) contain only one physical parameter, 
the friction Reynolds number 



lie = LUr/vQ , 



(3.5) 



and a set of geometry-dependent, dimensionless factors 
Gj , g n , and Sj . For the channel and pipe flows the factors 
Gj = VjjvQ and Sj are given in Table I. For simplicity in 
this paper we take g n = G n . 



B. Conservation laws and fluxes 

In the invisid (Tj = 7„ = 0), unforced (Vpj = 0) limit 
Eq. (3.3) conserves energy, £, linear momentum, V, and 
(quasi-) angular momentum, M.: 



71=1 



oo oc 

£ = aE^O^ + EK-l 2 )' (3 ' 6a) 

(3.6b) 
(3.6c) 



E^ 

3=1 



In general case, with non-zero Tj , jj , and Vpj the direct 
calculation of d£/dt, dV/dt, and dftA/dt with the help of 
Eq. (3.3) gives: 



d£ 
~dt 
dP 
dt 
dM 
~dT 



= e+-e- 



= m + - m 



(3.7a) 
(3.7b) 
(3.7c) 



Here e + , p + , and m + , are total influxes of respective 
integrals of motion, £, V, and A4, and e~ , p~, and m~, 
are their total rates of dissipations. Similar to Eqs. (3.6) 
these objects can be represented as Sj-weighted sums, 
J2j s j{- ■ ■ ) ; 01 the respective densities of partial influxes 
(e~j, pj~, and mi") and partial rates of dissipation in the 
j-zone: 



e+ = E ^ 4 - 4 = v ^ v l = v 'i - ( 3 - 8a ) 

oo 

P + = E'itf' P+ = V Pj = l, (3.8b) 

(3.8c) 



3=1 

OO 



m+ = ^ Sj m+ , m+ = , 

3=1 
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s = 



£j =F 1 \V J \ 2 + J2ln\u nj \ 2 , 

n=l 

(3.9a) 



i=i 

oo 



m = 



i=i 



(3.9b) 



(3.9c) 



The total influx of energy is exactly equal to the total 
linear momentum of the flow (due to our normalization, 
Vp =1) and the total influx of the linear momentum 
equals to unity: 

oo 

e+ = VpY^SjVj = VpV = V, (3.10a) 

i=i 

OO 

p+ = Vp^sj = X7p =1. (3.10b) 

i=i 

For the channel and pipe flows, the influx of the angular 
momentum, m + = 0. In the planar Couette flow, p + = 
and m + ^ 0. 

Notice, that the main part of the influx of the momen- 
tum is always flowing into the first few zones, where Sj 
have considerable values. For example, in the channel, 
~ 81% of the total p + flows into the first zone (occupy- 
ing si ~ 0.81 part of the cross-section area, see Table I), 
93% - into the first two zones and only about 1% to the 
fifth and all higher zones. As we show below, V- slightly 
decay with j, hence the influx of energy is even more 
confined to the first zones. 

For very high Reynolds numbers the dissipation of the 
conserved quantities occurs at large j (or n), ~ log 2 7?e. 
Therefore for log 2 Tie ^> 1 there exist an incrtial inter- 
val (generally speaking, different for the different quan- 
tities). Clearly, in the stationary case the influx of some 
conserved quantity have to be equal to its flux (in the 
shell space or via zones). For the energy, it is useful to 
analyze the energy flux in the shell space: 



€ S n £ 



= v. 



5 « = 2 ' 



(3.11) 



3=1 



where e n j is usual for the shell models (in the j-zone) 
expression for the flux of energy from nth to (n + 1) 
shell. 

The flux of the linear momentum, pj from j- to (j + 1)- 
zone is given by Eq. (2.37), which in the dimcnsionlcss 
form reads: pj — du 2 /2. Therefore, in the stationary 
case, for j in the inertial interval 



du 2 j /2 = p+ =p- = 1 



(3.12) 



To get the expression for the value of pj , valid for any j 
(not only in the inertial interval), we multiply Eq. (3.3a) 



by the j-zone area, Sj, and sum up from 1 to j. This 
gives the useful equation 



p j =du 2 j /2 = J2s i (l-T i V i ) , anyj, (3.13) 



that allows us to express the turbulent velocity Uj = Ujj , 
generating the turbulent energy cascade in the j-zone, 
via PM-velocities. 



C. Galilean invariance of the model 

The Galilean transformation to the reference system, 
moving with some velocity U$ in the streamwise direc- 
tion x, changes V(p, t) => G V(p, t) + x[/ . According to 
Eqs. (2.11a) and (2.28) it leads to the transformation: 



Vjit) ^ G Vj(t) + U 



(3.14) 



Clearly, the term Wj, Eq. (3.4c), describing the effect of 
the mean velocity on the turbulent subsystem, is Galilean 
invariant: 



V j+ i)u* 



3 ■ 



(3.15) 



In other words, the uniform velocity profile, V(p) = 
xxconst, does not affect the statistics of turbulence. This 
important invariance of the MZS-equations guarantees 
that for T^e — > oo the mean velocity profile in the core of 
the flow becomes uniform and the statistics of turbulence 
becomes space homogeneous, as expected. 



D. Asymptotic scale invariance and equations for 
the Turbulent Boundary Layer near flat plane 

Consider MZS-Eqs. (3.3) at very large Tie in the near- 
wall region, n, j > j* 3> 1 , in which the dimcnsionlcss 
parameters in Eq. (3.4) can be already replaced by their 
asymptotic values: 



G, 



9 j 9 , 
Oj = 2s j 



(3.16) 



This is definitely so for the channel and pipe flows, say, 
for = 3, see Tab. I. For large j the zone width is 
much smaller than the local curvature radius of the wall 
and the discussed situation corresponds to the case of 
Turbulent Boundary Layer (TBL) near the flat plane. 

For j > j* the dimcnsionlcss flux of linear momentum 
from j- to {j + l)-zonc, pj, is very close to 1 and much 
larger than the direct influx into j zone s.,p+ = Sj <C 1 
from the external pressure gradient. Therefore in this 
regime one can neglect in the RHS of Eq. (3.3a) Vp = 1 



16 



with respect to Wj and simplify Eq. (3.3) to the scale- 
invariant form: 



dVj 


Gk) 


dt 


-Re 


du n j 


9 K l 


dt 


-Re 


2 


dKjiVj 



^V,+dK 3 (u 2 3 _ 1 -U 2 j ). 



j'>3 

k = l/2s , n > j > j* , 



P + = rf<_ 2 /2 = rf<_ 1 /2 = l, 



(3.17a) 
(3.17b) 

(3.17c) 
(3.17d) 



where N n j has already the scale invariant form Eq. (3.4e). 
The "TBL boundary condition" at j = j*, (3.17d), pro- 
vides the influx of energy and mechanical momentum into 
TBL, Eq. (3.17). In initial Eq. (3.3) this role is played 
by the external pressure gradient Vp = 1. 

Equations (3.17) have additional [with respect to 
Eq. (3.3)] rescaling symmetry namely they remain un- 
changed under the transformation 



t 

-Re 



V 3 = Vj+30 ! 

u n j = u n +j a j+j a , 
t = t/2 ja , 
Ve = -Re/2 ja , 



(3.18) 



which corresponds (in the r-space) to the simultaneous 
rescaling of the outer scale L and the pressure gradient 
Vp in a way that leaves the value of the wall shear stress 
unchanged, \7pL — const. This symmetry of the equa- 
tions means that MZS-model describes the asymptotic 
universality of the near- wall turbulence for 7?e — > oo, and 
that the only relevant parameter in this regime is the to- 
tal influx of the momentum, p + , which is fixed by the 
boundary conditions (3.17d). 

In studies of the near-wall turbulence, it is customary 
to normalize the time and velocity units by the "inner" 
scales instead of the "outer" ones. Namely, we should 
use the viscous length-scale S = vq/U t = L/~Re instead 
of the outer scale L, and the corresponding timescale 
t$ = S/U T = r/Tvfe instead of r. As one can see, this 
rescaling corresponds to the choice j = log 2 7\!e in the 
transformation (3.18). In new units, the MZS-cquations 
have the same form as (3.17) with 7?e = 1: 



dVj_ 
dt 

du, 



-Gn 2 Vj + dK.j{u 2 _ Y — u 2 ) , 



dt 
A 



nj 2 



+- 



2 



dKjiVj - V j+1 )u* 



j'>3 



(3.19a) 
(3.19b) 

V-i' (N nf - N nj )] , 



where we omit tildes in the new variables. Equations 
(3.19) represent the MZS model for TBL near flat plane. 

Unlike Eqs. (3.17), in Eqs. (3.19) the zone and scale 
indices, j and n, can be both positive and negative: thus, 



Vq corresponds to the velocity of the scale <5, V-\ - to the 
structures of the scale 25, and so on. 

We can use general formulas (2.14a) to reconstruct the 
mean velocity profiles, V PM (p) . In this case the expres- 
sion for V PM (p) take on especially simple form, since all 
functions &j(p) already have a scale-invariant form. Tak- 
ing into account our rescaling to the "inner" units, we can 
write: 



V PM (y) = Rc 



,(2^^/26) 



(3.20) 



where y is the distance to the wall, < & U n(C) i s given by 
(2.20) and V PM (y) is measured in the units U T . 



IV. SOLUTION OF THE MZS EQUATIONS 
A. Laminar velocity profile and its instability 

1. Comparison of the full and PM-laminar profiles 

The simplest solution of MZS-cquations (3.3), in which 
all u n j — 0, corresponds to the laminar flow: 



-Re 



(4.1) 



Using Eq. (2.14a), one reconstructs the laminar profile of 
the PM-velocity 



(4.2) 



The laminar profile of the full velocity V°(p), satisfy the 
linear NSE, which in the dimcnsionless form (3.1) reads 



Ay°(p) + x% = o 



(4.3) 



According to Eq. (2.13) the PM-velocity is understood 
as the PM-projection of the full velocity: P PM {V(p)} = 
V PM (p). Therefore, we expect that the laminar PM- 
profilc (4.2) satisfy the equation, similar to (4.3): 



P PM {AV° M (p)}+ZKe 







(4.4) 



The proof of this equation is given in Appendix A. 

In the case P PM A = AP PM Eqs. (4.3) and (4.4) coin- 
cide, because P PM V PM = V PM . If so, the functions V°(p) 
and V pM (p) have to coincide due to uniqueness of solu- 
tions of Eq. (4.3) with zero boundary condition. How- 
ever, due to incompleteness of the PM-basis, P PM A ^ 
AP PM and hence the profiles V° (p) and V pM (p) are ex- 
pected to differ. In other words, Eq. (4.2) cannot recon- 
struct exactly the laminar profile V°(p). This is the price 
to pay for the incompleteness the PM-basis. Nevertheless 
the PM-basis is "full enough" to allow the reconstruction 
of any "physically possible" mean velocity profile with 
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The excitation of turbulence in, say, first zone will lead to 
a significant momentum flux from the 1st to the 2nd zone. 
As a result, the mean velocity Vi will increase, the 2nd 
zone velocity gradient Vi — V3 will increase too, and the 
instability condition for U2 will be therefore satisfied for 
smaller Reynolds numbers than predicted by Eq. (4.6). 
We will show later, that the real instability threshold for 
j-th zone is proportional to Tiej ex Kj, i.e. is much smaller 
than the "laminar" result (4.6). 



B. Wall-bounded flows in the approximation of 
near- wall eddies (a — b — c = 0). 



FIG. 5: Comparison of the exact (solid line) and PM- (dashed 
line) laminar profiles for the channel flow. The ratio of two 
profiles is shown by dotted line. Its values are marked on the 
right axis. 



good accuracy. As a first demonstration of this fact we 
compare in Fig. 5 the exact and the PM -laminar pro- 
files for channel flow, which differ only in few percents. 
This small loss of accuracy is insignificant with respect 
to a dramatic simplification of the calculation scheme 
for V(p, t): for large Tie the mean velocity field in PM- 
representation has N *=s \nTie significant coefficients Vj, 
while in the corresponding complete cell basis one has to 
account for ~ Tie N functions. 



2. Instability of the laminar flow at Tie — He c 



In the previous Sec. IV A we have found the laminar 
solution of the MZS-Eqs. (3.3) with u n j — and have 
shown that this solution becomes unstable at Tie = Tie cr 
with respect to excitation of the near wall velocity U\ = 
tin.. As the next step in this Section we analyze a much 
more general solution of Eqs. (3.3), that allows non-zero 
values of all near wall velocities, Ujj = Uj. The jj- 
eddies can be excited by a direct interaction with the 
PM- velocities, Vj. To prohibit the turbulent cascades, 
leading to excitation of other nj-eddies (with n > j), 
one can neglect in Eq. (3.3b) the interaction term J\f n j 
(i.e. put a = b = c = 0). In this case the MZS-Eqs. (3.3) 
take on a simple form: 



dt 

duj 

~~dT 



-TjVj +l + d Kj (v^ - u]) , (4.8a) 
-7^j + ^r(Vj-V j+1 )u*. (4.8b) 



The laminar solution (4.1) exists for all Tie. For large 
Tie, however, it becomes unstable with respect to excita- 
tion of turbulent near- wall eddies Uj. A simple analysis 
of the linearized Eq. (3.3b) (i.e. Af n j => 0) shows that 
the instability condition of j-eddy reads: 

7j <d(Vj-V j+1 )/2o-j . (4.5) 

Substituting here Vj from Eq. (4.1) one gets: 



9i K i 
Tie 



< 



d / 1 

2a.j VG, k 



r> ,-2 ) ^ ' 



'3 ^l™j 

which can be rewritten as the condition for Tie: 



Tie > Tie, 



9 3 Gj G j+ i 



Sj d{G J+ i - GjK]/K 2 j+1 ) 



It is clear that 7?e, oc k 3 / 2 , therefore, usually the first 

3 

unstable zone is the first one, j = 1. 

We want to stress, that Eq. (4.6) gives a correct insta- 
bility threshold only for the first unstable zone: 



(4.6) 



Tie cr = Tiei = k^ 2 . 



s 1 d(G 2 -G 1 



(4.7) 



In the stationary regime (d/dt = 0) only the finite num- 
ber (m — 1) > of turbulent velocities are non-zero, i.e., 
Uj = for j > m. This number depends on Tie. As it 
follows from Eq. (4.8a), the PM- velocities Vj for j > m 
coincide with the laminar ones: 

Vj = V° = TJ 1 - 2~ 2j Tie , j > m = m{Tie) . (4.9) 

In the unstable region j < to, as it is follows 
from Eq. (4.8b), the zone velocity difference, AVj = 
Vj — Vj + i coincide with it critical value, see Eq. (4.5): 



AVj = AV 3 C 



^m=Aj/-Re, (4.10) 



A 3 = -•l^'^j <> ■ 

This allows one to find all Vj via the last velocity in the 
unstable region, V m : 

Vj = V m + (B m -Bj)/Tie, (4.11) 
i-i 

Bj = 

i=l 
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Multiplying Eq. (4.8a) by Sj and summing up from j = 1 
to j = to, one finds the equation for V m with the solution: 



V m 

Drr. 



1 



0~m+l 



C n 



Tie 



D„ 



C m Tie 



i=i 
m 

i=i 



-) 2 m 



(4.12a) 
(4.12b) 

(4.12c) 



Note, that in spite of the large number of different pa- 
rameters (A m , . . . , D m ), all of them are provided with 
explicit expressions and can be easily evaluated. 

Now one finds the total linear mechanical momentum 
V (i.e. the total flux of the fluid) of the flow from its 
definition (3.6b) and Eqs. (4.9) and (4.12a) for Vj for 
j < to and for j > to: 



E 

j—m-\-l 



sV° 



E m -Re + F m /-Re. 



Era — 



F = 

± 771. 



(I -a. 



rn+l) 



E 



1 G A 

j—m-\-l J J 



(4.13a) 



(4.13b) 



m 

E^ 



Sj(B m - Bj)- 



(1 — a m+ i)D n 



So far we treated the last index to in the unstable zones 
as fixed. However, it depends on the Reynolds number 
Tie and is defined from the condition 

V m ~ Vm+1 < AV£ r = ^ , 

which ensures that the turbulent velocity u m will not be 
excited, u m = 0. This condition can be written as 



-Re 2 < 



A C 



1 



-Re m . (4.14) 



Actually it will be exactly (to — 1) excited turbulent ve- 
locities Uj, if 

Ve m -i <'Re<'Re m . (4.15) 

The analysis shows that this condition selects the index 
to, for which the total momentum V (4.13) has minimum 
value. 

Let us analyze the dependence (4.13a) of the total mo- 
mentum in the limit of extremely large Reynolds num- 
bers, ln7?e ^> 1, and, respectively, large critical zone 
numbers to. In this case we can use the scale-invariant 
limit for all parameters: 



En 



l E. 



2 m F, -Re m = 2 m R, (4.16) 



where E, F and R are some geometry-dependent con- 
stants. Now for Reynolds numbers inside the range (4.15) 




10 1 10 2 10 3 10 4 
Friction Reynolds number fie 



10 5 



FIG. 6: Total momentum of the flow vs. the Reynolds number 
Tie. Circles denote numerical data for the channel flow, d — 
0.046; solid lines - the analytical prediction (4.13) for different 
critical zone numbers m, m = 1, 2, . . . (from left to right). 



we can write in Eq. (4.13a) x{7te) = r Re/(2 m R), where 
1/2 < x(72e) < 1. Then V for large enough Tie [in the 
region of validity of the scale invariant limit (4.16)], is 
given by 



V(Ke) = x{-Re)ER 



x(Re)R 



(4.17) 



Since x(Tte) is bounded between 1/2 and 1, the flux (4.17) 
is bounded from above by one of the constants 



VCRe) < max | 



ER 



2F 
R 



ER 



R.i 



and does not grow infinitely for infinite Reynolds num- 
bers. In other words, neglecting the dissipation of energy 
in turbulent cascades, one concludes that with the given 
cross-section area of the flow and pressure gradient, the 
total flux of the fluid reaches some limit in spite of the 
infinite decrease of the kinematic viscosity. The numeri- 
cal and analytical calculations for the channel geometry, 
shown in Fig. 6, support this unexpected conclusion. 

The reason for this strange behavior is that the station- 
ary zone velocity differences, AV}, in the unstable region 
j < m are determined solely by the dissipation of turbu- 
lent eddies jj, Eq. (4.10), that go to zero as Tie — > oo (for 
finite j). Thus, in the limit Tie — > oo, AVj — > and Vj — > 
const. This conclusion is illustrated in Fig. 7. Clearly 
for large Tie first few velocities remain unchanged, and 
the only effect of increasing Tie is the shift of the "dissi- 
pative cutoff" towards the smaller scales. However, the 
total momentum is determined mainly by the first few 
zones and thus remains the same. 

This analysis shows importance of the turbulent cas- 
cade for the experimentally observed characteristics of 



19 





12-r 




1 u - 






w 


8- 


.9 




'o 




_o 


6- 






> 










4- 




2- 




0- 







\ \ \ \ \ 

WW R \ = 


- 

16000 


\ \ \ soooV 




\ \ 4000 \ \ 




\ 2000 \ \ \ 




1000 \ \ \ 









1 2 3 4 5 6 7 8 9 10 11 12 13 
Zone index j 



FIG. 7: Mean zone velocities Vj 
different Reynolds numbers Tie. 



the zone index j for 



real flows, i.e. infinite increase of the total flux, logarith- 
mic profiles, etc. The stationary condition AVj = AV" 
will hold even if one accounts for a turbulent cascade, 

T 

which adds some turbulent damping 7- for turbulent 
near- wall eddies instead of usual one jj. The turbulent 
damping does not vanish in the limit Tie — > 00. Thus, 
one obtains AVj — AVf 1 ^ 0, and the total momentum 
will infinitely increase, as it is shown in the following 
Sec. IV C. 



C. Wall bounded flows in the "turbulent-viscosity" 
approximation 

In this Section we show that even rough account for 
the turbulent cascade in the MZS-model already gives 
qualitatively correct analytical results. 



1. MZS equations in the turbulent viscosity approximation 

In the fully developed turbulent regime, the action of 
the interaction term M n j on the eddy be ap- 

proximately accounted for in the "turbulent viscosity" 
approximation in which the energy flux from the energy 
containing jj-eddies toward small scales in the j-zone is 
replaced by a nonlinear damping term, ensuring the same 
loss of their energy. Formally this can be done by replac- 
ing the nonlinear term Mjj in the full MZS-Eqs. (3.19) 

T 

by some effective turbulent damping term, 7 ■ : 



A/' 



« (a - c) . 



(4.18a) 
(4.18b) 



Here 7 • , is chosen as the turnover frequency of jj-eddies, 
with some dimensionless pref actor a. This pref ac- 
tor is evaluated in Eq. (4.18b) by equating the total rate 



of energy dissipation in the model system with the effec- 
tive turbulent damping (4.18a) and the energy flux to- 
ward small scales in the full shell model, see, e.q. Ref [3]. 

In the suggested effective turbulent damping approxi- 
mation, the MZS-Eqs. (3.19) take the form 



dV, 



-if - -TjVj + l + dKjiu 2 ^-^) 



dt 

duj 

~~dt 



(4.19a) 



— j- = -(j j + aK j \u j \)u j + — (V j -V j+ i)u 



2o-j 



(4.19b) 



These equations are different from the MZS-equations in 
the near-wall-eddy approximation, Eq. (4.8), by the only 
term oc a. 

For large Tie in the near wall region Eq. (4.19) have 
a more simple, scale-invariant form, in which all ge- 
ometry dependent factors are taken in their small scale 
limit (3.16): 



dt 

du 3 
dt 



Gk 2 

{9^1 1 ,\ dKj , Tr 

\liT ' " " ' ' ~ : l 



(4.20a) 
(4.20b) 



where Kj = V jls . 



2. Inertial interval solution 

In this section we consider the stable stationary solu- 
tion of Eq. (4.19). This solution is real. 

For log 2 Tie ^> 1, in the inertial interval of scales, all 
needed information contained in Eq. (4.19a) may be ob- 
tained from Eq. (3.13) for the momentum flux, pj, in 
which one can neglect the damping term cx 1 /Tie. This 
gives 



du) 



(4.21) 



where aj is defined by Eq. (3.4c). The substitution of 
Uj from Eq. (4.21) in Eq. (4.20b) gives an expression for 
Vj+i via Vj. This allows one to find Vj for any j outside 
the viscous region via some constant Vq: 



Vj = V Q 



— CTi+l 



d 3/2 

i=l 



(4.22) 



For j > j'*, (where is equal, say, 3) one can take in 
Eq. (4.22) Uj = 2sj <C 1. One concludes, that in the 
inertial interval the PM- velocities decrease linearly in j: 



Vj = V 



A - 

d 3/2 1^ 
i—1 



2V2 



d 3/2 



J > J* 



2- 



(Jiy/l — CTi+l 



(4.23a) 
(4.23b) 
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Here the geometry dependent constant A was found by 
comparison of Eqs. (4.22) and (4.23), and for different 
geometries it evaluates to: 



d 3/2 



2.7, 
2.4, 



channel , 
pipe . 



(4.24) 



The equation (4.23) can be directly obtained from the 
condition of the constancy of the momentum flux (3.12), 
that gives: 



const. 



(4.25) 



Then, in the stationary case, Eq. (4.20b) determines the 
j-independent difference Vj — Vj+i, necessary for keeping 
the amplitudes of the near-wall eddies at the constant 
level. This agrees with Eq. (4.23). 

In the derivation of Eqs. (4.22) and (4.23) we cancelled 
in Eq. (4.19b) Uj = u*, assuming that Uj ^ 0. At a given 
Tie this assumption is valid only for j < j . Here jo 
is the zone index of the last unstable zone, for which 
in Eq. (4.20b) d(lnuj)/dt > at Uj -> 0. To find j 
we consider Eq. (4.20b) as an equation for a continuous 
index q and set its RHS to zero with u q = 0: 



iHl - *(V 
Tie ~ 2 [Vq 



V q+ i) = 



V2a 



(4.26a) 



q(Tie) = log 2 7?e + log 2 2V2as/gVd .(4.26b) 

The index jo is then the integer part of q. 

The solution Eq. (4.22) has still an unknown constant, 
Vq. This constant can be found from the stationarity 
condition (3.12), 



The only 7?e-dcpendent factor here is the position of the 
viscous cutoff q(/Re), given by the Eq. (4.26b). Counting 
j from q(~Re), one has Tie independent PM- velocities. In 
the physical space this corresponds to the universality of 
the TBL profile measured in the "wall units" y + — yTZe. 

Notice, that in the derivation of Eq. (4.22), we ne- 
glected in Eq. (4.20b) the viscous damping term, gn^/TZe, 
with respect to the turbulent one, a\v,j\Kj. This approxi- 
mation fails near the viscous cutoff, because guj/Tie oc 2 J 
increases exponentially toward large j, while a\uj\ is ap- 
proximately constant. More detailed analysis shows, that 
Eq. (4.22) fails only for last two j before the cutoff and 
instead of sharp cutoff at j = jo there is a soft decrease 
of Vj in two zones near j = j , see Fig. 8. 

One can easily see, that the linear in j set of zone 
velocities (4.29) in the inertial interval corresponds to 
the log-profile of the mean velocity in the physical space: 
V(y) ~ \n(yTie/L). Actually, by a direct calculation one 
shows, that the log-profile 

V^(y) = —ln(^)+B, (4.30) 
corresponds to the zone velocities 



0-69 n „ „ 1.34 
{\og 2 Ke-j)+B 



(4.31) 



where 7 E ~ 0.58 is the Euler gamma constant. 

Thus, the MZS model describes the transition to the 
universal log-profile with the von-Karman constant 



£ 

3=1 



1 



(4.27) 



requiring that the total influx of the linear momentum, 
caused by the pressure gradient, must dissipate at the 
wall, due to the viscous friction. The product SjTj ~ 
2 J /Tie. Therefore, the sum in Eq. (4.27) is dominated by 
last few terms with j w j » j*. This allows one to use 
a more simple Eq. (4.23) instead of Eq. (4.22): 




2V2 



aj 



d 3/2 

2V2aj 

d 3/2 

2V2aq(Tie) 
d372 



(4.28a) 
(4.28b) 
. (4.28c) 



Together with Eq. (4.23a) this gives the dependence of 
Vj on j for the near wall region, that is geometry inde- 
pendent: 



2V2, 



d 3/2 



[q(Ve)-j] + 



Vd 



fJ 



V2aG 



j>3. (4.29) 



In 2 d 3 / 2 nn d 3 / 2 
k„ = — = ~ 0.25 



K 2V2 a 



a 



(4.32) 



D. Numerical analysis of the turbulent channel 
flow in the MSZ-model 

This Section is devoted to the numerical analysis of the 
MZS-modcl of the channel flow in the effective turbulent 
damping approximation, Eq. (4.19). For simplicity we 
adopt gj = Gj. PM-damping factor Gj was found ac- 
cording to Eq. (2.29) as the matrix elements of the vis- 
cous operator on the real part of the PM-basis &'j(y)- 
This gives 

Gj = 2 j+2 Sj . 

Parameters d — 4.6 x 10~ 2 and a = 5.9 x 10~ 3 were cho- 
sen to reproduce the experimental values of the universal 
constants k k w 0.4 and B w 5.2. 

In numerical analysis, 30 zones were sufficient to de- 
scribe the flows with Reynolds numbers up to fO 9 . In 
order to find the stationary solution of the full MSZ- 
Eqs. (4.19) we develop extremely stable and efficient it- 
eration procedure, based on the essential physics of the 



21 




10 15 20 
Zone index j 



o 



30 




-25 -20 -15 -10 -5 
"Wall" zone index j-\og 2 Re 



FIG. 8: PM-velocities Vj for different Reynolds numbers Tie as a function of the zone index j (left panel) and of the "near-wall" 
zone index (j ; — log 2 Tie) (right panel). The lines from lower to upper correspond to log 2 Tie = 9, 13, 17, 21, 25 and 29, respectively. 
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FIG. 9: Magnitude of the near- wall turbulent eddies \uj\ as 
a function of the zone index j for different Tie. The lines are 
marked as in Fig. 8. 



FIG. 10: Comparison of numerical solutions of the full MZS 
Eqs. (3.19) [black squares] and MZS Eqs. (4.19) in the approx- 
imation of turbulent viscosity (4.18) [empty circles]. Tie = 
10 6 . 



problem, see Appendix B. In spite of the huge Reynolds 
numbers, the accuracy better than 10 -6 was reached with 
about 100-200 iteration. Actually, using our approach 
(MZS-model and the iteration procedure) one can simu- 
late turbulent wall-bounded flows for arbitrary large Tie 
with very modest PC, even with XT486@40MHz, 8Mb 
RAM. 



1. Behavior of the zone velocities Vj and Uj 

In Fig. 8, left panel, we plot PM-velocities Vj for dif- 
ferent Tie from 500 to 5x 10 8 . One clearly sees the inertial 



interval, where Vj decrease linearly with j in agreement 
with Eq. (4.22). As we mentioned above there is a "soft" 
viscous cutoff, that involves last 2 zones. In order to 
demonstrate the phenomenon of universality, in the right 
panel of Fig. 8 we re-plot the same velocities Vj as func- 
tion of the "near- wall" zone index j — log 2 Tie as it is 
suggested by Eq. (4.29). There is a perfect collapse of all 
lines. Importantly, they collapse not only in the inertial 
interval, but also in the dissipative cutoff range. There is 
a non- negligible difference for the first two-three zones, 
j = 1,2, (see the blow up in the lower panel) that is 
caused by the non-negligible momentum influx in these 
zones. This deviation is described by Eq. (4.22). 
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Fig. 9 shows the magnitudes of the near-wall turbu- 
lent eddies, \uj\, for the same set of Tie's. Again, the 
numerical results are in full agreement with our theoret- 
ical understanding: the magnitudes of the near-wall ed- 
dies in the inertial interval are constant, as it required by 
Eq. (4.25) from the condition of the constancy of momen- 
tum flux. At the dissipative cutoff we observe a decrease 
of Uj involving 3-4 zones. A small decrease of Uj at the 
first few shells is caused by the non-zero momentum flux 
into these zones; in this region Uj with good accuracy 
can be found from 

pj =duj/2= 1-ctj+i . (4.33) 

To check how the approximation of turbulent viscos- 
ity, Eq. (4.18), affects the resulting mean velocity pro- 
file we compare in Fig. 10 the solutions of the full MZS 
Eqs. (3.19) (black squares) and MZS Eqs. (4.19) in the 
approximation of turbulent viscosity, Eq. (4.18) (empty 
circles). These solutions practically coincide except for 
three zones (j — 12, 13, 14) near the viscous cutoff. 
In these zones just a few shells are excited in the full 
model (3.19) and the observed difference is an artefact 
of discreetness of the scale space in the shell models (re- 
call, that in the shell model (3.19) the spacing parameter 
A = 2). In the case when the details of turbulent cascades 
near the wall are physically important (e.g., in turbulent 
flows with polymeric additives) one has to use more de- 
tailed shell representation of turbulent velocity field, say 
with A = V2 or even A = 2 1 / 4 , preserving the interaction 
range in the scale space unchanged. This modification 
of our model is under construction and will be published 
elsewhere. 

We expect that the more detailed models with A < 2 
will give the mean profile, closer to that, given by re- 
duced model with turbulent damping (4.19), in which 
the energy dissipation is affected by the scale spacing. In 
this paper we are not interested in details of turbulent 
cascades and will use further only Eqs. (4.19). 

2. Reconstruction of the mean velocity profile 

The equation (2.14a) reconstructs PM- velocity profile 
Vpm{v) from the set of Vj. For the case of the channel 
flow one has: 

oo 

v PM {y) = J2 v i*'j(y)> ( 4 - 34 ) 

3=1 
2 j -l 

$ i(y) = ^2 Pmhm-l{y), 

where p m and 4> m (y) for the channel geometry are given 
by Eq. (2.17). The results of such reconstruction for a 
set of Tie's are shown in Fig. 11. We see that recon- 
structed PM-profiles show all qualitative features of the 
"real" mean velocity profiles in the well-developed tur- 
bulent regime. Namely, the mean PM- velocity is almost 
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FIG. 11: Reconstructed PM -velocity profiles Vpm(v) for dif- 
ferent Reynolds numbers. 

constant in the main part of the flow [with centerline ve- 
locity, increasing with Tie as In(73e)]. As expected, all fall 
of V(y) occurs in a small region near the walls. It is also 
clear, that the width of this near- wall region decreases as 
Tie increases, in the same way as in real flows. 

Figure 11 shows, that the PM-profiles have some non- 
physical "wiggles" . They originate from the incomplete- 
ness of the PM-basis. Indeed, the PM-basis is con- 
structed by Eq. (2.4) from the complete (^>-basis with 
some prescribed m-dependence for each j function. As 
a result, the PM-expansion (4.34) can be understood 
as the Fourier expansion in 4> m {y) oc sm.{kim-xy) with 
discontinuities of the Fourier amplitudes at the "zone 
boundaries" m — 2 J , which produces the wiggles in the 
^-representation. This artifact of the model can be re- 
moved by different ways. The simplest one is to add some 
function V{y), orthogonal to all PM-basic functions (and 
thus having zero momentum), with the amplitude and 
the shape of which are determined from the problem of 
minimization of discontinuities in the spectrum. It can 
be shown that the result of such "smoothing" is fully 
acceptable for the most of purposes. 

Figure 12 displays the PM-velocities in the near-wall 
region y -C H in log- linear scale for different Tie's from 
~ 1.3 x 10 5 to ~ 5 x 10 s . The distance from the wall 
is measured in near- wall viscous lengths, y + = y/S = 
yRe/H. The collapse of profiles for different Reynolds 
numbers is evident. One can see the viscous sublayer 
(for y + < 10) and universal log-profile (for y + > 50). 



3. Comparison of the MZS universal mean velocity profile 
with experiment and DNS 

In Fig. 13 we compared the reconstructed MZS uni- 
versal mean velocity with the DNS results in chan- 
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FIG. 12: Collapse of the mean- velocity profiles in the near- 
wall units. Different line (from below to above) correspond 
to Tie = 2 17 ~ 1.3 x 10 5 ,2 21 ,2 25 and 2 29 ~ 5 x 10 8 
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FIG. 13: Comparison of the reconstructed MZS universal 
mean velocity profile V(y + )/U r , (solid line in the region 
1 < V + < 10 3 ) with the DNS results in channel, available 
in Ref. [15] for y + < 10 2 (dashed line) and with the mea- 
surements in pipe, taken from Ref. [16] for y + > 10, (empty 
circles) . 

nel, Ref. [15] and with the laboratory measurements at 
Tie up to 3.5 x 10 7 , presented in Ref. [16]. The DNS data 
are available for y + < 100, while the experimental data 
are obtained for y + > 10. As expected, in the overlap 
region, 10 < y + < 100 both results collapse. As we ex- 
plained, our MZS model reproduce the asymptotical log- 
profile (4.30), and the parameters of the model, a and d 
were chosen to give known values of k k and B, that pa- 
rameterize Eq. (4.30). Therefore, as expected, the MZS 
profile, displayed in Fig. 13 as solid line, coincide with 
the experimental date in the region of large enough y + . 



The point is that the MZS dependence V PM {y + ) practi- 
cally coincides with the DNS and experimental data in 
all region of y + . It means that the MZS model correctly 
describes the basic physics that affect the mean veloc- 
ity profile in the universal near wall region of turbulent 
boundary layer near the flat plane. The MZS descrip- 
tion of the viscous and the buffer layers does not require 
adjustable parameters. 

4- Reconstruction of the profile of the energy dissipation 

The total dissipation rate ej in j-zone, Eq. (3.9a), can 
be splited into the dissipation rate in the PM-velocity 
subsystem ej and the dissipation rate in the turbulent 
subsystem i j : 

ej = ej+ej, (4.35) 

i<j 

In the turbulent viscosity approximation the effective 
damping is given by 

7ofj =7j +««jKI • (4-36) 

Figure 14 displays as empty circles the dissipation den- 
sity in turbulent subsystem, normalized by the near-wall 
lengthscale, i.e. ej /Tie as a function of the "wall zone 
number" defined here as 

j+ =j-log 2 7& + 2 . (4.37) 

The black squares in this figure show the PM-velocities 
in the "wall-zone" representation, i.e. values of Vj+ for 
the same Tie = 5 x 10 s . The solid line is the log- plot 
of the reconstructed (from this set of Vj) PM-velocity in 
the physical space, i.e. V PM (y + ) vs. log 2 y + , shown from 
above. As one sees, the solid line goes very close to the 
black squares, as it should be according to the interpreta- 
tion of the PM-velocities Vj+ as a "physical" velocity at 
some point within j + -zone, that explained in Sec. II B 6. 
As it is clear from Fig. 14, the wall zone numbers (4.37) 
are chosen to give a very simple "correspondence rule" : 

Wall zone index j + <^ Wall distance y + = 2~ j+ . (4.38) 

This approach can be used to restore of spatial distribu- 
tion of various turbulent characteristics. In particular, 
we can understand e~ + , presented in Fig. 14, as e~(y + ) 

with y + = (2~- 7+ ). As expected, the dissipation rate 
(normalized, as in Fig. 14, by Tie) is 7?e-independent as 
Tie — * oo. In this sense, the results in Fig. 14 can be con- 
sidered as universal. As expected, at large lie the main 
energy dissipation occurs in a narrow near-wall region, 
y + < 40. The dissipation rate at the wall in our model 
is equal to 

e | wail = Hm ej/Tie « 0.21 , Tie > 10 4 , (4.39) 
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FIG. 14: Mean velocities Vj+ (black squares) and turbulent 
dissipation rates (empty circles) in the wall-zone represen- 
tation j + , introduced by Eq. (4.37). Solid line - reconstruc- 
tion of the PM- velocity profile (in the wall units) V PM (y + ) 
vs. distance y + in the wall units, shown above. Tie = 5 x 10 s . 

which is reasonably close to the result for smaller Tie, 
available in the DNS of Ref. [17]: 

e| wall = 0.166, 7&- 200. (4.40) 

Notice, the difference in Tie's and that our result (4.39) 
is obtained in the simple model with just two adjustable 
parameters, chosen to adjust very different characteris- 
tics of the flow, the mean velocity profile. This allows one 
to consider the reasonable correspondence of Eqs. (4.39) 
and (4.40) as an argument in favor of our simple MZS 
model. 



5. Tie-dependence of the global flow characteristics 

Clearly, our approach allows one to evaluate various 
global characteristics of the turbulent wall bounded flows. 

The first example is the 7?e-dependence of the total 
momentum (i.e. total flux) of the flow, shown in Fig. 15. 
There is a laminar regime for Tie < 17, and a developed 
turbulent regime for, say Tie > 100. The "period-2" os- 
cillations are an artifact of the model, caused by the dis- 
cretization with the spacing parameter A = 2. These 
oscillations, however, are small and should be ignored. 
In principle, they can be removed in more "advanced" 
versions of the MZS models with few variables in each 
j-zone, Vj, a , Uj !<7 , responsible for "cr-sub-zones of the j- 
zone. 

The second example is the 7\!e-dependence of the total 
energy of the system, as well as parts of the energy, con- 
taining in the mean and turbulent subsystems, shown in 
Fig. 16. 

Last, but not least, example is the 7?e-dependence of 
the total energy dissipation and that for the V- and u- 
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FIG. 15: Dependence of the total momentum of the flow V 
from the Reynolds number Tie in the channel flow. 
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FIG. 16: Dependence of the total energy (A), the energy of 
the mean V-subsystem (B), and the energy of the turbulent 
ii-subsystem (C) from the Reynolds number Tie. 

subsystem, shown in Fig. 17. Total energy dissipation 
is equal to the total energy influx, i.e. to the total mo- 
mentum of the flow. However, the distribution of energy 
dissipation between two subsystems is very interesting. 

One clearly sees that some flow characteristics, like the 
energy of turbulent subsystem and the energy dissipation 
in the mean flow subsystem remain finite for vanishining 
fluid viscosity, vq. In the same time, other characteris- 
tics, such as the total linear momentum, the energy of the 
mean flow sub-system and the rate of energy dissipation 
increase infinitely [like ln(l/Vo)], i.e. demonstrate a phe- 
nomenon of viscous anomaly. The MZS model clearly 
demonstrate that the physical reason for that anomaly 
is the separation in the physical space of the external 
forcing and the friction: the external pressure gradient, 
that accelerate the flow, acts on the whole cross-section 
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FIG. 17: Dependence of the total energy dissipation (A), en- 
ergy dissipation in the mean V^-subsystem (B), and energy dis- 
sipation in the turbulent u-subsystem (C) from the Reynolds 
number Tie. 

area of the flow, while the friction force, that prevents 
the mean velocity from infinite growth, acts only on the 
walls. To be able to maintain the constant flux of the 
linear mechanical momentum toward the wall, the am- 
plitudes of the near- wall eddies of all scales, Uj, must 
be j-independent. This immediately leads to a linear in- 
crease of Vj with j, decreasing from the viscous cutoff 
value, jo, toward the beginning of the cascade, j = 1. 
Thus the value of V\ is proportional to the total number 
of the cascade steps, jo — 1- Decrease of uq to a half of 
its value adds one more step in inertial interval of the 
momentum cascade. This leads to an increase of V\ by 
~ In 2. In the dimensional units this corresponds to the 
increase V\ by ~ V T In 2 ~ ln2yTVp . 



V. SUMMARY 

We develop a multi-zone shell (MZS) model for the 
wall-bounded turbulent flows in the piecewise approxima- 
tion, dividing cross-section area into set of N ~ log 2 Tie 
j-zones. In each zone the turbulence is assumed to be ho- 
mogeneous and is described in the framework of a shell 
model equation for the "turbulent" shell velocities u n j(t). 
The mean velocity is described by an additional set of N 
zone variables Vj (t) , that allow to reconstruct the mean 
velocity profile with the help of a specially designed PM- 
basis. 

The MZS model 

• conserves the actual integrals of motion of the origi- 
nal NSE of the problem, energy, linear and angular mo- 
menta; 

• respects Galilean and "asymptotic" scale invariance, 
NSE type of nonlinearity. 

• in the relatively simple and analytically transparent 



manner describes the basic physical phenomena in the 
wall bounded flows for a huge interval of Reynolds num- 
bers. They include: 

1. the laminar velocity profile for TZe < TZe cr ; 

2. its instability at TZe = TZe cr , 

3. intermediate, non-universal mean velocity profile at 
moderate Tie; 

4. universal profile for TZe ^> TZe CI in the viscous sub- 
layer, buffer layer and log-law region; 

5. spatial distribution of turbulent activity, of the rate 
of energy dissipations, etc. 

• allows additional adaptation of the MZS equations for 
first few (energy containing) shells to the particular flow 
geometries (like the channel, pipe, Couette flows, etc.) 
that may be based on the stability analysis of the lam- 
inar regime or some other specific geometry determined 
information. This should improve description of the flow 
for moderate Tie < 1000. 

• may be generalized on the case of viscoelastic turbu- 
lent flows (by adding additional shell variables for the 
polymeric additives), on the case of particle laden sus- 
pensions, etc. 
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APPENDIX A: EQUATION FOR THE LAMINAR 
PM -PROFILE 

Here we present a proof that the laminar PAf-profile, 
V° M (P): given by Eq. (4.2), satisfies Eq. (4.4). 

To this end, we recall that the damping parameters 
GjK? are defined as matrix elements: 

Gi«? = - a 7 1 (*j» A *i) • ( A1 ) 

Substituting this definition into (4.1), (4.2), we obtain 

KM) = -*E7^y- < A2 > 

Then, the Laplacian of is equal to 

j \ 3> 3/ 
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Acting on this Laplacian with the PM-projector (2.12), 
we obtain 

where one must distinguish the Reynolds number, Tie, 
from the notation of the real part of something, Re{. . . }. 
Note that by construction, different PM-functions 3>j be- 
long to the different subspaces of the eigen-functions of 
the Laplace operator. Also, the real and imaginary 
parts of these functions are orthogonal, and thus 

(* fc ,A*;.) = A fcj (*;.,A*;.) . (A5) 

As a result, all non-diagonal terms in (A4) vanish, and 
we obtain 

P PM {AO = (A6) 

3 

Finally, using property (2.11a), we obtain 

P PM {&V°J = -'ReSt, (A7) 
which actually is Eq. (4.4). 

APPENDIX B: ITERATION PROCEDURE FOR 
SOLVING MZS EQUATIONS IN THE 
EFFECTIVE TURBULENT DAMPING 
APPROXIMATION 

Consider the MZS-Eqs.(4.19) for some large Tie. As 
we discussed, there are two regions of j, namely j < 
jmax and j > j max . In the "turbulent region", j < j max , 
the laminar solution, Uj = 0, is unstable with respect to 
excitation of turbulent amplitudes Uj and thus, in the 
stationary regime, Uj ^ 0. For d > 0, without loss of 
generality all Uj can be taken real and positive definite. 
In the turbulent region the stationary velocities Vj and 
Uj satisfy the equation: 

TjVj = l + dnjiu^-u*), (Bla) 

jj + aKjUj = jL(Vj-V j+1 ), (Bib) 

uj > , j < j max . 

In the "laminar region", j > j max : 

Vj = 1/Tj , Uj =0, j > j max . (B2) 



One can substitute Vj from Eq. (Bla) to Eq. (Bib) 
and get an equation, connecting the triad Uj—i, Uj, and 
Uj+i. Unfortunately, a direct iteration procedure in this 
"straightforward" equation is unstable and do not con- 
verge to the stationary solution. 

To find the stable, stationary solution of Eqs. (Bl) nu- 
mericcaly, we develop a stable iteration procedure, which 
is based on the physical understanding of these equations 
as describing the momentum flux from j = 1 toward large 
j in the "interaction triads" Vj, V}+i and Uj. Denote as 
vj p \ Vj^ and u^ the solution of following three alge- 
braic equations on the j-step in the p-th iteration run: 

TjV^ = l + d K ,{[^M-f] 2 }, 

r.+iv^i = i + d K3+1 {[uff-[u^f}, 

li+otKjuf = A-[V^-V}1\}, (B3) 
uf > , for j max > j > 1 , 4 P) = ■ 

As the "initial condition" at p = we take the (unstable) 
laminar solution: 

vf> = r-, nf) = o. 



In the first step, j = 1, of each iteration run one takes 
u ( v ] = 0. Finding v} p) , vj p) , and u (p) one takes in 
Eq. (B3) j = 2, finds V^ p) , V^ p) and u { f ] and so on until 
on some jo step one gets negative (or complex) solution 
foru^. It means that this amplitude is stable and has 
to be taken zero, — 0. Accordingly, j = j max . For 
all j > j max one takes the laminar solution: V^ = 1/Tj, 
uf = 0. 

After that one begin the next, p + 1, iteration run, 
starting again from its first step, j = 1. It can be shown, 
that the velocities V^ form monotonically decreasing 
sequence with increasing p and are always positive. Since 
a limited from below, monotonically decreasing sequence 
always have some finite limit, this proves the convergence 
and stability of our iteration scheme. The calculations 
show, that for Reynolds numbers Tie < 10 9 velocities Vj 
and Uj converge (with accuracy about 10~ 6 ) after 100 — 
200 iterations runs. 
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